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Abstract 


The  aim  of  the  current  research  has  been  to  contribute  to  the  “aerodynamic  forcing 
function”  aspect  of  the  high  cycle  fatigue  problem.  In  this  regard  density  wakes  were  iden¬ 
tified  as  a  potential  new  source  of  turbomachinery  blade  vibration.  In  order  to  characterize 
the  density  wake  induced  force  and  moment  fluctuations  a  two-dimensional  computational 
study  was  conducted  to  simulate  the  passage  of  density  wakes  through  a  cascade  blade 
row.  Both  inviscid  incompressible  flows  and  viscous  compressible  flows  ( M oo  =  0.15  to 
Moo  =  0.87  and  Reynolds  number  Re(c,  U0 o)  ~  700,000)  were  considered. 

The  inviscid  flow  simulations  show  a  pair  of  counterrotating  vortices  in  the  blade  passage 
and  a  flux  of  density  wake  fluid  to  the  blade  suction/pressure  surface  as  low/high  density 
wakes  convect  through  the  blade  row.  This  fluid  flux  is  found  to  be  the  key  mechanism 
which  affects  the  blade  static  pressure  distribution  and  hence  the  blade  force  and  moment 
response. 

The  viscous  flow  simulations  have  uncovered  additional  mechanisms  for  blade  force  and 
moment  fluctuations:  (1)  periodic  vortex  shedding  at  the  blade  trailing  edge,  (2)  convecting 
separation  bubbles  on  the  blade  suction  surface  (due  to  the  density  wake  -  blade  boundary 
layer  interaction)  and  (3)  axial  deflection  of  the  blade  passage  shock  wave  position  (due  to 
the  density  wake  -  shock  wave  interaction). 

Simple  functional  relationships  have  been  derived  to  quantify  the  force  and  moment 
fluctuations  for  each  flow  simulation.  A  simple  cascade  flow  model  has  also  been  developed 
to  determine  parametric  trends  in  the  blade  force  and  moment  fluctuations  with  varying 
density  wake  properties  and  compressor  geometries. 
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1  Background 

Increased  operational  requirements  and  increased  thrust  to  weight  ratios  have  led  to  higher 
mean  and  fluctuating  stresses  in  components  of  modern  turbomachinery.  This  has  increased 
the  likelihood  of  encountering  high  cycle  fatigue  (HCF)  failure  in  fan,  compressor  and  tur¬ 
bine  blades. 

The  large  number  of  parameters  and  the  wide  occurrence  of  HCF  producing  conditions 
over  the  engine  operating  regime  imply  that  structural  integrity  must  be  evaluated  in  an 
extremely  large  number  of  situations.  An  additional  implication  is  that  it  is  difficult  to 
extract  general  guidelines  for  high  cycle  fatigue  prevention  because  of  the  high  dimension¬ 
ality  of  the  parameter  space  that  must  be  explored  [1].  Further  complexity  is  introduced 
by  the  diversity  of  local  phenomena,  e.g.  tip  leakage  flows,  unsteady  shock  motion  and 
local  separation  that  are  characteristic  of  turbomachinery  flows.  The  lack  of  HCF  design 
principles  which  consider  these  flow  phenomena  is  exemplified  in  the  following  statement 
taken  from  a  recent  HCF  workshop  held  at  the  MIT  Gas  Turbine  Laboratory  [1]:  “...forced 
blade  response  is  not  currently  predictable,  and  structural  design  and  analysis  for  high  cy¬ 
cle  fatigue  situations  have  not  advanced  beyond  the  early  concepts  of  the  fatigue  limit,  the 
Goodman  diagram  and  Miner’s  rule. 

The  aim  of  the  current  research  is  to  focus  on  aerodynamic  characterization  of  blade 
force  response  to  help  increase  the  predictive  capability  of  HCF  failure1 .  Several  forced  vi¬ 
bration  “sources”  exist  in  turbomachinery.  In  particular  viscous  wakes  from  upstream  blade 
rows  and  potential  flow  effects  due  to  rotor-stator  interaction  have  received  a  lot  of  atten¬ 
tion  (Kemp  and  Sears  [13],  Kerrebrock  and  Mikolajczak  [14],  Manwaring  and  Wisler  [17], 
Valkov  [27])  2.  The  forced  vibration  source  considered  in  this  study  are  convecting  density 
wakes.  Density  wakes  can  enter  the  engine  from  ground  ingested  hot  air,  steam  ingestion 
during  carrier  launches  and  exhaust  gas  ingestion  from  forward  firing  weapons  [21].  The 
difference  in  temperature  between  the  blades  and  the  surrounding  fluid  can  also  generate 
density  gradients  particularly  in  the  downstream  wakes  of  blades.  Incomplete  and  non- 
uniform  combustion  can  also  introduce  density  gradients  into  the  flow  entering  turbines. 


2  Technical  Approach 

Computational  fluid  dynamics  (CFD)  is  the  “tool”  used  in  this  research  to  investigate  the 
forces  and  moments  induced  by  convecting  density  wakes.  CFD  is  a  relatively  inexpensive 
and  convenient  method  to  investigate  flow  phenomena  compared  to  experimental  inves¬ 
tigations  in  wind  tunnels.  In  particular  the  flow  field  properties  at  any  location  in  the 
computational  domain  can  be  conveniently  determined  and  an  overall  “picture”  of  the  flow 

^ote  that  forced  vibrations  are  distinct  from  self  excited  vibrations  due  to  blade  flutter  and  dynamic 
stall. 

2Note  these  studies  are  primarily  concerned  with  compressor  performance  and  the  corresponding  effects 
on  unsteady  blade  loading  have  not  been  considered. 
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can  be  generated  to  help  identify  specific  flow  features.  The  flow  geometry  and  free  stream 
conditions  are  also  easily  changed  within  a  CFD  simulation  compared  to  an  experimental 
facility  where  arbitrary  changes  in  flow  geometry  are  generally  not  feasible. 

The  use  of  CFD  is  constrained  by  available  computational  resources  however.  This 
research  has  therefore  been  limited  to  two-dimensional  unsteady  flows  with  a  single  com¬ 
pressor  blade  row  geometry.  While  turbine  blades  are  subjected  to  larger  density  non¬ 
uniformities  (due  to  hot-streaks  from  the  combustor  and  from  blade  cooling),  compressor 
blades  are  considered  here  since  they  are  more  susceptible  to  HCF  failure. 

The  compressor  blade  profile  used  in  this  study  is  the  General  Electric,  Low  Speed 
Research  Compressor  (LSRC)  Stator-B  blade.  The  geometry  of  this  blade  is  described  in 
Section  14.3. 

The  density  wakes  considered  convect  along  the  axial  direction  and  have  density  gradi¬ 
ents  directed  only  in  the  axial  direction.  Discussion  is  focused  on  low  density  wakes  (wake 
densities  lower  than  free  stream  density)  which  are  more  common  in  compressor  blade  pas¬ 
sages.  To  help  isolate  individual  flow  features  a  single  density  wake  is  convected  through  the 
blade  passage  in  all  simulations.  The  density  variation  inside  the  density  wake  is  specified 
to  be  sinusoidal.  This  variation  is  considered  to  be  a  representative  case. 

To  investigate  the  wide  parameter  space  of  density  wake  -  cascade  blade  row  interac¬ 
tions  a  simple  cascade  model  is  also  developed.  The  model  uses  a  combination  of  singularity 
solutions  and  a  proportional  constant  determined  from  the  inviscid  CFD  results. 

Subsequent  sections  in  this  report  are  organized  as  follows  : 

•  Physical  description  of  the  origin  of  density  wake  induced  blade  forces. 

•  Theoretical  background  which  includes  the  derivation  of  non-dimensional  scaling  re¬ 
lationships. 

•  Results  from  the  inviscid  flow  simulations 

•  Results  from  the  viscous  flow  simulations. 

•  Results  from  the  cascade  flow  model. 

•  Appendices: 

—  Description  of  the  viscous  flow  solver,  computational  grid  and  non-dimensionalization. 

-  Blade  force  and  moment  fluctuation  profiles. 

1.  Baseline  response  and  boundary  layer  profiles. 

2.  Primary  and  secondary  response  -  time  histories  and  functional  relationships. 

3.  Primary  response  compressibility  scaling. 

-  Cascade  model  source  code. 
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Figure  1:  Density  wake  convecting  through  a  compressor  blade  row. 


3  Physical  Origin  of  Unsteadiness 

An  analogy  to  the  passage  of  a  density  wake  through  a  compressor  blade  row  can  be  found 
in  the  atmosphere;  low  density  (high  temperature)  air  rises  to  higher  altitudes  where  the 
pressure  is  lower  and  remains  there  because  of  force  equilibrium.  Similarly,  high  density 
(low  temperature)  air  sinks  to  regions  of  high  pressure  nearer  to  the  earth. 

Now  consider  the  passage  of  a  density  wake  through  a  cascade  blade  row  as  illustrated 
in  Figure  1.  Assume  the  wake  has  a  lower  density  than  the  free  stream  density.  As  the 
wake  moves  through  the  cascade,  the  low  density  fluid  migrates  towards  the  suction  side 
(low  pressure  side)  of  the  blades  by  the  action  of  centrifugal  forces.  To  satisfy  mass  conser¬ 
vation  the  surrounding  higher  density  fluid  is  subsequently  displaced  toward  the  pressure 
side  of  the  blade.  This  relative  motion  of  low  and  high  density  fluids  generates  a  pair  of 
counterrotating  vortices  in  the  blade  passage.  The  low  density  fluid  directed  toward  the 
blade  suction  surface  and  the  associated  counterrotating  vortices  convect  through  the  blade 
passage  together  with  the  density  wake. 

The  blade  pressure  distribution  is  influenced  by  the  impact  of  the  low  density  fluid  on 
the  blade  surface.  The  blade  force  and  moment  coefficients  therefore  change  with  time 
during  passage  of  the  density  wake. 


4  Theoretical  Background:  Marble’s  Linearized  Analysis 

For  the  case  of  a  density  wake  convecting  through  a  compressor  blade  row  as  shown  in 
Figure  1  the  density  gradients  and  pressure  gradients  are  misaligned  by  almost  90  degrees. 
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[a] 


Figure  2:  Lift  coefficient  fluctuation  during  passage  of  a  density  discontinuity  over  a  flat  plate.  A 
is  the  position  of  the  density  discontinuity  as  it  convects  along  the  flat  plate.  The  flat 
plate  lies  between  |A|  <  1. 


This  misalignment  allows  for  vorticity  production  in  the  flow3.  For  a  continuous  density 
distribution  in  two-dimensional,  inviscid,  incompressible  flow  this  vorticity  satisfies  the 
linearized  relation, 


( d  d\  1 

(at  +  “si ) “  =  x  v"  <*) 

If  the  density  gradient  (Vp)  is  large  (zeroth  order),  the  convected  vorticity  u>  is  of  the  same 
order  as  the  pressure  field  [18]. 

The  first  study  of  density  gradients  as  a  source  of  flow  unsteadiness  was  conducted  by 
Marble  [18].  He  performed  a  linearized  potential  flow  analysis  for  a  flat  plate  at  an  angle 
of  attack  a  encountering  a  density  discontinuity.  If  the  fluid  is  treated  as  incompressible 
and  the  velocity  disturbances  caused  by  the  airfoil  are  small  compared  to  the  free  stream 
velocity,  the  density  field  can  be  expressed  as  p(x  -  ut,  y)  [18].  The  results  of  Marble’s  anal¬ 
ysis  are  shown  in  Figure  2  and  Figure  3.  In  these  Figures  A  is  the  position  of  the  density 
discontinuity  as  it  convects  over  the  flat  plate.  The  flat  plate  lies  between  |A|  <  1. 


3The  generation  of  vorticity  due  to  the  misalignment  of  density  gradients  and  pressure  gradients  is  often 
referred  to  as  “baroclinic  torque”. 


9 


Figure  3:  Moment  coefficient  fluctuation  during  passage  of  a  density  discontinuity  over  a  flat 
plate.  A  is  the  position  of  density  discontinuity  as  it  convects  along  the  flat  plate.  The 
flat  plate  lies  between  |A|  <  1. 


The  dotted  lines  correspond  to  the  quasi-steady  results  while  the  solid  lines  correspond 
to  the  unsteady  results.  Initially  the  effect  of  the  density  discontinuity  is  to  reduce  the 
local  lift.  This  is  a  consequence  of  a  downwash  velocity  field  which  precedes  the  arrival  of 
the  discontinuity.  This  is  followed  by  a  rapid  rise  in  lift  as  the  discontinuity  convects  over 
the  leading  edge.  This  is  caused  by  an  upwash  velocity  field  behind  the  discontinuity.  A 
gradual  relaxation  of  the  perturbation  occurs  as  the  density  discontinuity  convects  further 
downstream.  The  final  steady  lift  scales  with  the  ratio  of  the  density  across  the  discontinu¬ 
ity.  The  moment  coefficient  shown  in  Figure  3  also  reflects  these  events  in  local  loading. 

Marble’s  linearized  analysis  provides  a  basic  understanding  of  the  parameters  involved 
in  this  problem.  In  particular  the  density  parameter  p*, 


P2  -  Pi 

p  —  - 

P2  +  Pi 


(2) 


is  shown  to  be  a  key  parameter  in  the  unsteady  loading.  Re-writing  Equation  1  using  the 
density  parameter  p*  and  the  non-dimensionalized  vorticity  Cb  =  uc/Uoo  gives, 


Du 

~Df 


(3) 
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where, 


- 

p 

p  — 

(P2  -  Pl)Ul 

r  = 

2s, 

c 

p 

p  — 

Pi  +  P2 

vw  = 

wV 

Vfc  = 

hV 

Equation  3  suggests  the  non-dimensional  wake  width  w/c  and  non-dimensional  blade  spac¬ 
ing  h/c  to  be  additional  key  parameters. 


11 


x/c 


Figure  4:  Perturbation  velocity  vectors  during  passage  of  a  density  wake  of  width  0.2c  and  density 
ratio  0.5  through  the  NACA4F  blade  row.  The  flow  is  inviscid  and  incompressible. 
r  =  0.65. 

5  Results  From  The  Inviscid  Flow  Simulations 

The  inviscid  flow  simulations  capture  the  essential  flow  physics  of  the  density  wake  -  blade 
row  interaction  described  in  Section  3.  In  particular  the  passage  of  low  density  wake  fluid 
to  the  blade  suction  surface  and  the  formation  of  counterrotating  vortex  pairs  within  the 
blade  passage  are  observed.  The  perturbation  velocity  vectors  in  Figure  4  obtained  during 
passage  of  a  density  wake  of  width  0.2c  and  density  ratio  0.5  clearly  illustrate  these  flow 
features. 

The  transport  of  the  low  density  wake  fluid  to  the  blade  suction  surface  constitutes  a 
fluid  flux  which  can  be  expressed  mathematically  as, 

rw/c 

fluid  flux  —  f  pvjdx  (4) 

Jo 

where  p  is  the  wake  density,  vj  is  the  velocity  of  the  wake  fluid  and  dx  is  a  increment 
in  length  parallel  to  the  axial  flow  direction,  p  and  vj  are  functions  of  x.  vj  and  non- 
dimensional  wake  width  w/c  are  also  functions  of  time  r. 

The  impingement  of  the  low  density  fluid  on  the  blade  suction  surface  results  in  a  local 
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increase  in  the  blade  static  pressure.  As  the  density  wake  convects  downstream  the  static 
pressure  increase  moves  further  aft  along  the  blade  suction  surface.  This  effect  can  be  seen 
in  the  suction  surface  pressure  distribution  plotted  in  Figure  5.  A  corresponding  decrease 
in  static  pressure  occurs  at  the  blade  pressure  surface  due  to  the  “suction”  effect  of  the 
counterrotating  vortices.  The  net  result  is  a  local  reduction  in  the  static  pressure  difference 
across  the  blade  and  a  time  varying  fluctuation  in  the  blade  force  and  moment. 

The  magnitude  of  the  force  and  moment  fluctuations  increases  as  the  fluid  flux  to  the 
blade  surface  increases.  Hence  larger  wake  widths  and  smaller  density  ratios  increase  the 
magnitude  of  the  blade  force  and  moment  fluctuations  4.  Figure  6  shows  the  time  vary¬ 
ing  azimuthal  force  coefficient  and  moment  coefficient  during  passage  of  a  range  of  density 
wakes.  Initially  a  reduction  in  the  azimuthal  force  occurs  as  the  density  wake  convects  over 
the  front  half  of  the  blade.  This  is  followed  by  a  gradual  return  to  steady  state  as  the 
density  wake  leaves  the  blade  trailing  edge.  Similarly  an  increase  in  the  counterclockwise 
moment  occurs  as  the  density  wake  convects  over  the  front  half  of  the  blade.  As  the  density 
wake  passes  the  blade  mid  chord,  the  moment  decreases  back  towards  steady  state.  The 
shape  of  these  force  and  moment  profiles  are  roughly  common  over  the  range  of  density 
wake  widths,  0.1  <  w/c  <  0.4,  and  density  ratios  0.25  <  P2/ Pi  <  2.0  considered. 

The  amplitude  of  the  maximum  fluctuations  in  the  force  and  moment  coefficients  were 
found  to  scale  linearly  for  small  wake  widths  w/c  and  small  density  parameter  p*  (see 
Figure  7  and  Figure  8).  This  scaling  is  given  by, 


Cymax  Cyr] 
Cy-mean 


CfTlmax  Cmr 

C  TTlmean 


(5) 

(6) 


4  Smaller  density  ratios  are  associated  with  larger  density  gradients  which  generate  stronger  vorticity  and 
higher  v/. 
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surface  distance  /  chord 


Figure  5:  Suction  surface  pressure  distribution  during  passage  of  a  density  wake  width  0.2c  and 
density  ratio  0.5. 


-0  5  0  0  5  1  1-5  2 

Time 


Figure  6:  Fluctuation  in  (a)  azimuthal  force  coefficient  and  (b)  moment  coefficient  (positive 
counter-clockwise  about  the  mid-chord)  during  passage  of  density  wakes.  A Cy  = 
Cymax  —  Cymcan  A  Cm  —  Ct/Ijua j  “*  Cnimenn  W  here  Cy  mean  —  0.75,  Cm  mean  —  -0.13. 
Reproduced  from  Ramer 


Figure  7:  Maximum  fluctuation  in  the  azimuthal  force  coefficient  (from  steady  state)  during 
passage  of  density  wakes  through  the  NACA4F  blade  row.  A Cy  —  ( Cymax  — 
C y mean) / Cy mean  *  Reproduced  from  Ramer. 


Figure  8:  Maximum  fluctuation  in  moment  coefficient  (from  steady  state)  during  passage  of  den¬ 
sity  wakes  through  the  NACA4F  blade  row.  A  Cy  =  (Cymax  —  Cy  mean)  /  Cy  mean- 
Reproduced  from  Ramer. 
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5.1  Conclusions 

•  The  passage  of  density  wakes  through  a  compressor  blade  row  results  in  a  time  varying 
fluctuation  in  the  blade  force  and  moment  coefficients.  This  fluctuation  is  a  result  of 
the  localized  reduction  of  static  pressure  difference  across  the  blade  due  to  a  flux  of 
density  wake  fluid  to  the  blade  surface.  During  passage  of  low  density  wakes  the 
wake  fluid  is  directed  toward  the  blade  suction  surface.  This  reduces  the  blade  force 
coefficient  and  increases  the  counter-clockwise  moment  coefficient.  The  opposite  is 
true  during  passage  of  a  high  density  wake.  In  this  case  the  wake  fluid  is  directed 
toward  the  blade  pressure  surface. 

•  Parametric  studies  show  the  amplitude  of  the  maximum  fluctuation  in  blade  force  and 
moment  coefficients  to  have  the  following  functional  relationship: 

A Cy  =  f{w/c, p*,  Cy{mean )) 

A  Cm  =  f(w/c,p*,Cm(mean)) 

The  effect  of  blade  spacing  h/c  is  included  in  Cy[mean )  and  Cm(mean). 

•  The  shape  of  the  force  and  moment  coefficient  fluctuations  are  common  over  the 
range  of  density  wake  widths  0.1  <  w/c  <  0.4  and  density  ratios  0.25  <  P2/P1  <  2.0 
considered. 

•  For  w/c  <  0.2  and  —0.2  <  p*  <  0.2,  A  Cy  and  A  Cm  scale  linearly.  Increasing 
non-linearity  is  observed  for  larger  wake  widths  and  density  ratios. 
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Convective  Time  Scale  (r) 

Location  of  density  wake 

0.0 

wake  l.e.  intercepts  blade  l.e. 

0.1 

wake  t.e.  intercepts  blade  l.e. 

1.3 

wake  t.e.  passes  blade  t.e. 

Table  1:  The  location  of  a  density  wake  at  different  times  during  passage  through  the  LSRC 
cascade  blade  row.  w/c  —  0.2,  P2/P1  =  0.50,  Moo  =  0.15. 


6  Results  From  The  Viscous  Flow  Simulations 

The  viscous  flow  simulations  were  conducted  at  Reynolds  number  Reilly,  c)  ^  700,000 
and  Mach  numbers  ranging  from  =  0.15  to  =  0.87.  The  purpose  of  these  sim¬ 
ulations  was  to  investigate  the  impact  of  blade  boundary  layers  and  blade  passage  shock 
waves  on  the  density  wake  blade  row  interaction.  The  important  features  of  the  viscous 
flow  simulations  are  best  described  using  the  time  varying  blade  force  and  moment  profiles 
shown  in  Figure  9.  The  force  and  moment  profiles  can  be  divided  into  3  different  regions 
depending  on  the  position  of  the  density  wake  within  the  blade  passage  (see  Table  1).  The 
characteristic  flow  features  in  each  response  region  are  summarized  below. 


1.  Baseline  response  region.  This  response  corresponds  to  all  times  when  the  density 
wake  is  upstream  of  the  compressor  blade  leading  edge.  The  force  and  moment  fluc¬ 
tuations  observed  here  are  a  consequence  of  vortex  shedding5. 

2.  Primary  response  region.  This  response  corresponds  to  when  the  density  wake  is 
within  the  cascade  blade  passage.  The  characteristic  flow  features  here  are  the  fluid 
flux  directed  to  the  blade  surfaces  and  the  associated  counterrotating  vortices6.  For 
flow  at  Mqo  =  0.87  blade  passage  shock  wave  fluctuations  are  also  observed. 

3.  Secondary  response  region.  This  response  corresponds  to  all  times  after  the  density 
wake  leaves  the  cascade  blade  trailing  edge.  A  separation  bubble  is  observed  on  the 
blade  suction  surface  during  this  response  7 . 

6.1  Vortex  Shedding  From  Cascade  Blade  Rows 

The  baseline  response  vortex  shedding  was  observed  for  Reynolds  numbers  in  the  range 
Re(UoQ,c)  =  600,000  —  800,000.  A  summary  of  the  force  and  moment  fluctuation  am¬ 
plitudes  and  frequencies  are  included  in  [26].  Experimental  results  from  Lawaczeck  and 
Heinemann  [16],  Boldman  et.al.  [3]  and  Parker  [19]  have  indicated  vortex  shedding  at  sim¬ 
ilar  flow  conditions  and  from  similar  cascade  geometries. 

5 Note  vortex  shedding  is  a  result  of  the  viscous  flow  environment  and  is  not  a  consequence  of  the  passage 
of  density  wakes . 

6Note  these  are  the  same  features  as  for  inviscid  incompressible  flows. 

7A  qualitative  description  of  this  flow  feature  has  been  included  in  the  September  1998  Progress  Re¬ 
port  [26]. 
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Baseline  Response  Primary  Response  Secondary  Response 


Figure  9:  Fluctuations  in  (1)  azimuthal  force  coefficient,  (2)  axial  force  coefficient  and  (3)  moment 
coefficient  (positive  clockwise  about  the  mid  chord)  during  passage  of  a  density  wake 
of  width  0.2c  and  density  ratio  P2/P1  =  0.5.  3  distinct  regions  can  be  identified  in  the 
response.  Moo  =0.15. 

6.2  Primary  Response  Force  And  Moment  Fluctuations 

The  time  varying  blade  force  and  moment  fluctuations  for  the  viscous  flow  simulations  are 
plotted  in  Figure  54  through  Figure  105  in  Section  16.  The  magnitude  of  the  force  and 
moment  fluctuations  are  plotted  separately  in  Figure  10  through  Figure  13. 

The  change  in  magnitude  the  force  and  moment  fluctuations  with  density  wake  width 
and  density  ratio  (at  a  given  flow  Mach  number)  follow  the  same  trends  as  the  inviscid  flow 
predictions,  i.e.  the  magnitude  increases  with  increase  in  both  wake  width  and  wake  density 
ratio.  This  is  as  expected  as  the  mechanism  which  consists  of  the  fluid  flux  to  the  blade 
surface  is  responsible  for  the  change  in  blade  static  pressure  distribution  for  both  viscous 
and  inviscid  flow  simulations. 

In  the  viscous  flow  simulations  the  vortex  shedding  induced  force  and  moment  fluctu¬ 
ations  are  superimposed  on  the  density  wake  induced  force  and  moment  fluctuations.  In 
most  cases  the  maximum  fluctuation  in  the  primary  response  can  be  easily  distinguished 
from  the  baseline  response.  However  as  the  free  stream  Mach  number  increases  the  baseline 
fluctuations  increase  in  magnitude  and  can  sometimes  exceed  the  primary  response  fluc¬ 
tuations.  This  is  particularly  true  during  the  passage  of  density  wakes  with  large  density 
ratios  (P2/P1  =  0.75)  at  high  free  stream  Mach  numbers.  A  clear  distinction  between  the 
baseline  response  and  the  primary  response  cannot  be  seen  for  these  cases.  The  maximum 
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Figure  10:  Maximum  fluctuation  in  the  force  and  moment  coefficients  as  functions  of  density 
wake  width  and  density  parameter  p* .  The  straight  lines  joining  data  points  are  to 
help  aid  clarity.  Moo  =  0.15. 


fluctuation  in  the  primary  response  is  nevertheless  obtained  by  identifying  the  maximum 
(or  minimum)  point  on  the  response  curves  corresponding  to  the  time  when  the  density 
wake  is  contained  within  the  blade  passage. 


6.2.1  Compressibility  Scaling 

Comparison  of  the  results  in  Figure  10  through  Figure  13  show  that  in  general  the  force 
and  moment  fluctuations  for  a  given  density  wake  width  and  density  ratio  to  increase 
with  freestream  Mach  number.  This  increase  can  be  scaled  using  the  Prandtl-Glauert 
compressibility  factor  as  follows, 


Cy0(max) 

-  Cy{max)  x  y/l  - 

(7) 

Cx0(max) 

=  Cx(max)  x  yj  1  — 

(8) 

Cm0(max) 

—  Cm(max)  x  yjl  — 

(9) 
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Figure  11:  Maximum  fluctuation  in  the  force  and  moment  coefficients  as  functions  of  density 
wake  width  and  density  parameter  p* .  The  straight  lines  joining  data  points  are  to 
help  aid  clarity.  Moo  =  0.53. 


Figure  12:  Maximum  fluctuation  in  the  force  and  moment  coefficients  as  functions  of  density 
wake  width  and  density  parameter  pm .  The  straight  lines  joining  data  points  are  to 
help  aid  clarity.  A/oc  =  0.63. 
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Figure  13:  Maximum  fluctuation  in  the  force  and  moment  coefficients  as  functions  of  density 
wake  width  and  density  parameter  p* .  The  straight  lines  joining  data  points  are  to 
help  aid  clarity.  Moo  —  0.87. 


where  the  “o”  subscript  indicates  the  incompressible  value.  Figure  106  through  Figure  117 
in  Section  17  compares  the  force  and  moment  coefficients  scaled  in  this  manner  at  each 
Mach  number.  The  following  conclusions  can  be  made  from  these  plots. 

•  For  wake  widths  of  0.1c  the  maximum  fluctuation  in  the  azimuthal  force  coefficient 
for  all  compressible  tests  compares  to  within  5%  of  the  corresponding  incompressible 
inviscid  result.  For  wake  widths  of  0.2c  the  deviations  from  the  incompressible  result 
increase  from  5%  to  a  maximum  of  20%  for  |p*|  >  1/3.  Similarly  for  wake  widths  of 
0.4c  deviations  of  up  to  40%  can  be  seen  for  |p*|  >  1/3. 

•  For  wake  widths  of  0.1c  the  maximum  fluctuation  in  the  axial  force  coefficient  for 
Mqo  =  0.53,  Mqo  —  0.63  and  M^  —  0.87  deviates  up  to  a  maximum  of  35%  from 
the  corresponding  result  at  M^  —  0.15  for  all  p*.  For  wake  widths  of  0.2c  and  0.4c  a 
deviation  of  more  than  50%  can  be  seen  between  the  compressible  results. 

•  For  wake  widths  of  0.1c  and  0.2c  the  maximum  fluctuation  in  the  moment  coefficient 
deviates  up  to  35%  between  the  M <*>  =  0.15  and  the  incompressible  calculation.  For 
wake  widths  of  0.1c  there  is  a  maximum  deviation  of  60%  for  p*  <  0  and  a  maximum 
deviation  of  100%  for  p*  >  0.  For  larger  wake  widths  the  deviations  exceed  100%. 

Clearly  the  Prandtl-Glauert  factor  is  not  adequate  to  describe  the  compressibility  effects 
on  the  maximum  fluctuation  in  the  force  and  moment  coefficients.  Deviations  from  this 
factor  increase  with  increasing  Mach  number,  increasing  wake  width  and  increasing  density 
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ratio.  Reasonable  scaling  for  the  azimuthal  force  coefficient  can  be  achieved  however  for 
small  wake  widths  (w/c  —  0.1). 

6.2.2  Shock  Wave  Fluctuations 

The  viscous  flow  simulation  at  M0 0  =  0.87  is  characterized  by  a  shock  wave  0.25c  down¬ 
stream  of  the  blade  leading  edge  as  shown  in  Figure  14. 

The  passage  of  a  density  wake  through  the  blade  row  causes  fluctuations  in  the  posi¬ 
tion  of  the  shock  front.  This  in  turn  causes  changes  in  the  static  pressure  rise  across  the 
shock  wave.  The  maximum  deflection  of  the  shock  wave  as  a  function  of  the  density  wake 
width  and  density  parameter  is  plotted  in  Figure  15.  During  passage  of  low  density  wakes 
(p*  <  0.0)  the  shock  front  is  initially  displaced  upstream  towards  the  blade  leading  edge.  As 
the  density  wake  moves  further  aft  the  shock  wave  returns  back  towards  the  unperturbed 
position.  Note  for  a  high  density  wake  (wake  density  higher  than  free  stream  density)  the 
shock  wave  moves  downstream.  In  general  the  shock  wave  displacement  increases  with  den¬ 
sity  wake  width  and  density  parameter  \p*\-  The  maximum  deflection  of  the  shock  wave  is 
found  to  vary  almost  linearly  for  values  of  density  parameter  p*  <  —0.15.  The  data  points 
for  density  parameter  p*  =  —0.60  and  wake  widths  0.2c  and  0.4c  indicate  the  last  recorded 
position  of  the  shock  wave  before  it  is  suppressed  by  the  higher  temperature  (lower  Mach 
number)  within  the  density  wake  [26]. 

The  effect  of  shock  motion  on  (1)  the  static  pressure  rise  across  the  shock  and  (2)  the 
blade  azimuthal  force  coefficient  is  illustrated  in  Figure  16.  Notice  the  maximum  reduc¬ 
tion  in  azimuthal  force  coefficient  corresponds  to  the  point  of  maximum  static  pressure  rise 
across  the  shock. 

While  there  is  additional  contribution  to  changes  in  the  blade  static  pressure  distribu¬ 
tion  from  the  unsteady  shock  wave  motion,  this  contribution  is  such  that  the  blade  force 
and  moment  fluctuation  magnitudes  increase  with  density  wake  width  and  density  ratio. 
Hence  there  is  no  change  in  the  trends  for  the  force  and  moment  fluctuations  between  flows 
with  shock  waves  and  without  shock  waves. 


6.3  Secondary  Response  Force  And  Moment  Fluctuations 

The  magnitude  of  the  secondary  response  fluctuations  depend  on  the  fluctuation  of  the 
suction  surface  flow  separation  point.  This  fluctuation  in  separation  point  is  due  to  the 
formation  and  subsequent  downstream  convection  of  the  separation  bubble.  The  maximum 
upstream  motion  of  the  separation  point  is  plotted  in  Figure  17  through  Figure  20  for 
varying  density  wake  widths  and  density  ratios  for  each  Mach  number  flow.  In  general  the 
separation  point  moves  further  upstream  as  the  wake  density  ratio  decreases  (for  a  high 
density  wake,  P2/P1  >  1>  the  separation  point  moves  further  upstream  as  the  density  ra¬ 
tio  increases).  An  increase  in  density  wake  width  also  drives  the  separation  point  further 


22 


Figure  16:  The  change  in  (1)  static  pressure  difference  across  the  shock  wave,  (2)  azimuthal  force 
coefficient  and  (3)  blade  shock  location  during  passage  of  a  density  wake  width  0.1c 
and  density  ratio  0.25.  Moo  =  0.87. 


upstream  however  there  is  an  upper  limit  ( w/c  —  0.2)  beyond  which  further  increase  in 
density  wake  width  does  not  further  effect  the  separation  point  motion.  This  is  particularly 
evident  in  the  M0 0  =  0.15,  =  0.53  and  =  0.63  results  for  density  ratios  p2/Pi  <  1- 

The  magnitude  the  secondary  response  force  and  moment  fluctuations  are  plotted  in 
Figure  21  through  Figure  24  for  varying  density  wake  widths  and  density  ratios.  The 
following  conclusions  can  be  made  from  these  plots. 

•  The  magnitude  of  the  force  and  moment  fluctuations  follow  the  same  trends  as  the 
separation  point  fluctuation  with  changes  in  density  wake  width  and  density  ratio. 
Note  how  the  curves  for  wake  widths  of  0.2c  and  0.4c  almost  collapse  together.  This 
corresponds  to  the  same  trend  observed  for  the  separation  point  fluctuation  in  Fig¬ 
ure  17  through  Figure  20.  This  shows  a  clear  link  between  the  separation  point  motion 
and  the  magnitude  of  the  secondary  response  fluctuations. 

•  The  magnitude  of  the  secondary  fluctuations  normalized  by  the  baseline  fluctuation 
magnitudes  decrease  with  increasing  free  stream  Mach  number.  This  however  is  a 
consequence  of  the  increased  magnitude  of  baseline  fluctuations  with  increased  Mach 
number.  In  general  the  secondary  fluctuations  for  all  coefficients  without  normaliza¬ 
tion  increases  with  Mach  number. 

The  formation  of  a  separation  bubble  on  the  blade  suction  surface  temporarily  suppresses 
the  vortex  shedding  at  the  blade  trailing  edge.  This  results  in  a  lack  of  high  frequency  con- 
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Figure  17:  The  maximum  change  in  the  suction  surface  separation  point  from  the  mean  baseline 
position  as  a  function  of  wake  width  and  density  ratio.  Moo  =0.15 


Figure  18:  The  maximum  change  in  the  suction  surface  separation  point  from  the  mean  baseline 
position  as  a  function  of  wake  width  and  density  ratio.  Moo  =  0.53 
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Figure  19:  The  maximum  change  in  the  suction  surface  separation  point  from  the  mean  baseline 
position  as  a  function  of  wake  width  and  density  ratio.  Moo  —  0.63. 


Figure  20:  The  maximum  change  in  the  suction  surface  separation  point  from  the  mean  baseline 
position  as  a  function  of  wake  width  and  density  ratio.  A/oo  =  0.87. 
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Figure  24:  The  maximum  fluctuation  in  the  blade  azimuthal  force,  axial  force  and  moment  coef¬ 
ficients  in  the  secondary  response  region  for  varying  density  wake  widths  and  density 
ratios.  A/oc  =  0.87. 
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Figure  25:  Discrete  fourier  transform  of  the  baseline  force  and  moment  coefficients.  dft(X)  is 
the  discrete  fourier  transform  of  the  time  signal  X.  p  =  non-dimensional  frequency. 
Moc  =  0.87. 


tent  in  the  force  and  moment  fluctuations.  This  is  evident  by  comparing  the  frequency 
spectrum  plots  obtained  separately  for  the  baseline  response  region  (Figure  25)  and  the 
secondary  response  region  (Figure  26).  Visual  examination  of  the  time  histories  of  the  force 
and  moment  fluctuations  for  all  the  parametric  tests  (Figure  54  through  Figure  105  in  Sec¬ 
tion  16)  also  indicates  this  result. 


6.4  Conclusions 

•  The  cascade  blade  force  and  moment  fluctuations  observed  in  viscous  compressible 
flows  can  be  divided  into  3  main  regions  based  on  the  position  of  the  density  wake: 

1.  The  baseline  response.  This  corresponds  to  when  the  density  wake  is  upstream 
of  the  blade  leading  edge. 

2.  The  primary  response.  Response  when  the  density  wake  is  contained  within  the 
blade  passage. 

3.  The  secondary  response.  Response  after  the  density  wake  leaves  the  blade  trailing 
edge. 

•  The  baseline  response  is  a  consequence  of  vortex  shedding  at  the  blade  trailing.  The 
frequency  of  blade  force  and  moment  fluctuations  are  found  to  correspond  with  the 
vortex  shedding  frequency  and  the  magnitude  is  found  to  increase  with  flow  Mach 
number  [26]. 
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Figure  26:  Discrete  fourier  transform  of  the  secondary  response  force  and  moment  coefficient 
fluctuations,  w/c  =  0.2,  p2/p\  =  0.5,  Moo  =  0.87.  dft(X)  is  the  discrete  fourier 
transform  of  the  time  signal  A',  p  =  non-dimensional  frequency. 


•  The  magnitude  of  the  primary  response  force  and  moment  fluctuations  scale  with  the 
magnitude  of  the  fluid  flux  directed  to  the  blade  surfaces  (similar  to  the  inviscid  flow 
force  and  moment  fluctuations).  At  M0 0  =  0.87  shock  wave  deflection  is  also  found  to 
contribute  towards  the  force  and  moment  fluctuations.  The  magnitude  of  the  fluid  flux 
and  the  deflection  of  the  shock  wave  both  scale  with  (1)  the  non-dimensional  density 
wake  width  w/c  and  (2)  the  density  parameter  p* .  The  trends  for  the  maximum 
fluctuation  in  the  force  and  moment  coefficients  in  viscous  flows  (primary  response) 
are  therefore  identical  to  the  trends  observed  in  the  inviscid  flows. 

•  The  magnitude  of  the  primary  response  force  and  moment  fluctuations  increase  with 

flow  Mach  number.  For  small  density  wake  widths  ( w/c  =  0.1)  the  magnitude  of  the 
azimuthal  force  coefficient  fluctuations  arc  found  to  scale  with  the  Prandtl-Glaucrt 
compressibility  correction  factor  \J\  —  . 

•  The  magnitude  of  the  force  and  moment  fluctuations  in  the  secondary  response  are 
typically  2-10  times  larger  than  the  vortex  shedding  induced  baseline  fluctuations. 
These  additional  fluctuations  are  caused  by  a  separation  bubble  on  the  blade  suction 
surface  formed  by  the  density  wake  -  boundary  layer  interaction.  The  magnitude  of 
these  fluctuations  scale  with  the  maximum  fluctuation  in  the  blade  suction  surface 
separation  point  which  in  turn  scales  with  (1)  the  normalized  density  wake  width  w/c 
and  (2)  the  density  ratio  pi/ P\- 

•  The  separation  bubble  is  found  to  temporarily  suppress  vortex  shedding  for  2-3  con- 
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vective  time  scale  units  after  the  density  wake  has  left  the  blade  trailing  edge.  The 
vortex  shedding  frequencies  are  consequently  absent  in  the  frequency  response  of  the 
secondary  fluctuations. 
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7  Cascade  Flow  Model 


CFD  flow  solvers  are  useful  for  investigating  specific  fluid  flow  phenomena  however  they  are 
less  useful  for  conducting  parametric  studies.  This  is  due  to  the  excessive  computational 
run-times  involved4 5 * * 8.  Furthermore  investigation  of  changes  in  flow  geometry  is  tedious  as 
a  new  computational  grid  must  be  generated  for  each  simulation.  A  simple  flow  model 
which  incorporates  the  essential  fluid  dynamic  features  can  however  greatly  reduce  the  time 
required  to  determine  parametric  trends  and  can  also  help  assess  the  relative  importance 
of  flow  geometry  variables  on  a  specific  output  variable. 

A  cascade  flow  model  was  therefore  developed  in  this  study  in  parallel  with  the  CFD 
simulations  to  help  determine  the  density  wake  induced  blade  force  and  moment  fluctua¬ 
tions  in  an  economical  manner.  The  results  from  the  viscous  compressible  CFD  simulations 
in  particular  indicate  a  variety  of  fluid  phenomena  associated  with  the  passage  of  density 
wakes;  these  cannot  be  easily  modeled  at  the  present  time.  As  a  first  step  we  develop  a  sim¬ 
ple  model  aimed  at  characterizing  trends  in  the  force  and  moment  fluctuations  for  inviscid 
incompressible  flows  only.  The  assumptions  and  equations  used  in  the  model  are  presented 
below  together  with  the  results  from  a  parametric  analysis.  The  source  code  for  the  model 
is  listed  in  Section  18. 


7.1  Modeling  Assumptions 

The  model  is  focused  on  the  density  wake  induced  force  and  moment  fluctuations  in  inviscid 
incompressible  background  flows.  A  combination  of  potential  flow  singularity  solutions  and 
a  constant  of  proportionality  determined  from  the  inviscid  CFD  results  is  used  to  develop 
the  model.  In  addition  the  model  uses  “quasi-steady  aerodynamic”  assumptions  and  is 
therefore  at  best  valid  for  low  reduced  frequencies  (/i  <  1.0  [4]).  The  following  additional 
assumptions  are  made  to  further  simplify  the  model: 

1.  The  flow  is  two  dimensional. 

2.  The  cascade  blades  are  unstalled  so  the  flow  always  follows  the  blade  surface. 

3.  Effects  due  to  blade  camber  and  thickness  are  neglected  so  that  the  cascade  blades 
can  be  represented  as  flat  plates. 

4.  The  density  wake  can  be  represented  as  a  row  of  counterrotating  vortices  which  convect 
with  the  mean  flow. 

5.  The  influence  of  shed  vorticity  can  be  neglected  in  the  transient  response. 

A  schematic  of  the  model  is  shown  in  Figure  27.  The  basic  approach  involves  solving 

for  the  bound  vorticity  on  the  cascade  flat  plates  such  that  the  flow  remains  tangent  to  the 

8The  investigation  of  time-accurate  flow  phenomena  in  particular  requires  extensive  computational 
resources. 
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Figure  27:  Model  schematic  indicating  flat  plate  cascade  and  counterrotating  vortex  pairs. 


plate  surface.  The  flat  plates  are  divided  into  N  discrete  panels  each  with  a  bound  vortex 
at  the  1/4  chord  location.  There  are  no  changes  in  the  bound  vortex  circulation  strength 
on  similar  panels  on  different  plates  (i.e.  zero  interblade  phase  angle).  The  flow  tangency 
condition  is  applied  at  the  3/4  chord  location  of  each  panel.  The  cascade  flat  plates  are 
inclined  at  stagger  angle  (3\  and  have  space-chord  ratio  a.  The  free  stream  velocity  is  in¬ 
clined  at  an  angle  a  to  the  flat  plates. 

The  convecting  density  wake  is  modeled  by  a  row  of  counterrotating  vortex  pairs.  This 
approach  is  motivated  from  the  inviscid  CFD  flow  visualization  results  which  indicate  the 
formation  of  counterrotating  vortex  pairs  as  the  density  wake  interacts  with  the  cascade 
pressure  field  (see  Figure  4).  The  counterrotating  vortices  are  staggered  at  angle  (3\  (equal 
to  the  blade  stagger  angle)  and  are  convected  parallel  to  the  flat  plates  at  a  fixed  vertical 
offset  h.  The  density  wake  width  is  represented  by  the  spacing  d  of  the  counterrotating 
vortices.  The  circulation  Ta  and  Tb  associated  with  each  of  the  vortices  in  the  counterro¬ 
tating  vortex  pair  is  prescribed  as  a  function  of  space  and  time. 


7.2  Induced  velocities 

As  stated  above  the  basic  approach  used  in  the  model  involves  solving  for  the  circulation 
strength  of  each  flat  plate  vortex  panel  such  that  the  flow  tangency  condition  is  satisfied 
on  the  blade  surface.  The  first  step  in  this  approach  requires  the  evaluation  of  the  induced 
velocity  field  of  (1)  the  bound  vortices  and  (2)  the  convecting  counterrotating  vortex  pairs. 
For  the  case  of  an  isolated  flat  plate  the  velocity  w (airfoil)  induced  by  a  single  bound  vortex 
of  circulation  strength  T  at  a  distance  r  from  the  vortex  is  given  by, 

r 

w(airfoil)  =  - —  (10) 

2nr 

For  the  case  of  a  cascade  of  flat  plates  with  bound  vorticity  at  N  discrete  vortex  panels  the 
induced  velocities  from  a  row  of  vortices  extending  from  —  oo  to  ~f  oo  must  be  summed.  The 
velocity  w (cascade)  induced  at  any  location  C  by  a  row  of  vortices  strength  T*  located  at 
the  same  relative  position  (  ona  discrete  vortex  panel  is  given  by 

w (cascade)  =  ^  el®'coth  — -  (11) 

where  j3\  is  the  plate  stagger  angle  and  T  is  the  plate  spacing.  This  expression  is  derived  in 
Kaufmann  [12].  The  real  and  imaginary  parts  of  Equation  11  yield  the  induced  velocity  in 
directions  parallel  to  and  normal  to  the  flat  plate.  For  the  counterrotating  vortices  which 
convect  a  distance  h  above  the  flat  plate,  an  equivalent  expression  can  be  obtained  for  the 
induced  velocities.  The  velocity  field  w (vortices)  induced  by  the  front  row  of  vortices  in  the 
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counterrotating  pair  (Figure  27)  each  with  circulation  strength  and  stagger  angle  (3\  is 
given  by, 


w  (vortices)  = 


-TfiT el0'coth  |  (tt^-^  +  t7r/l)  } 


(12) 


7.2.1  Circulation  Strength  of  The  Counterrotating  Vortices 

The  circulation  strength  of  the  vortices  in  the  counterrotating  pair  is  assumed  to  be  a 
function  of  the  vorticity  oj  (of  the  vortex),  an  effective  area  and  time  r  ,  i.e. 


r  =  f{u,Aeff,T ) 


(13) 


The  vorticity  u  is  assumed  to  change  from  zero  to  a  fixed  magnitude  as  the  vortices  convect 
through  the  cascade.  This  change  in  vorticity  is  governed  by  the  interaction  of  the  wake 
density  gradient  and  the  blade  pressure  gradient  according  to  Equation  1.  A  hyperbolic 
tangent  functional  approximation  is  assumed  to  represent  this  change  in  vorticity.  By 
interchanging  the  time  dependency  r  with  the  vortex  position  £  (by  assuming  a  constant 
convection  speed)  the  change  in  counterrotating  vortex  circulation  strengths  T  can  now  be 
expressed  as, 


T  =  0.5  x  rmax  x  (1  +  tanh(s +  52)) 


(14) 


where  Tmax  is  a  constant  dependent  on  the  density  wake  properties  and  the  cascade  pres¬ 
sure  distribution.  The  variables  S\  and  S2  are  chosen  such  that  the  tank  function  varies 
between  —3  and  +3  as  the  counterrotating  vortex  location  £  varies  between  predetermined 
end  points  l\  and  Z2.  These  end  points  are  obtained  by  considering  the  upstream  and  down¬ 
stream  decay  of  the  cascade  pressure  field  as  a  function  of  the  blade  row  stagger  angle, 
space-chord  ratio  and  blade  leading  edge  radius.  For  values  of  f  less  than  l\  the  circulation 
strengths  are  set  to  zero  while  for  values  greater  than  I2  the  circulation  strengths  are  set  to 
r max-  The  —3  to  +3  range  of  the  tank  function  is  chosen  to  ensure  a  smooth  transition  of 
TA  and  F#  to  the  fixed  end  values  (see  Figure  28). 

rmaI  in  Equation  14  is  assumed  to  contain  the  influence  of  the  blade  row  pressure 
gradient  and  the  density  gradient  together  with  an  effective  area  as  follows, 


Tmax  —  / ^  Pi  Ae //) 


(15) 
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Figure  28:  Change  in  circulation  strength  of  the  convecting  vortices.  v=  si£  *f  $2  where  £  is  the 
location  of  the  convecting  vortex  and  s\  and  «2  are  constants  which  ensure  the  tank 
function  is  evaluated  between  —3  and  +3. 


In  the  simplest  case  the  functional  relation  for  rmax,  Vp,  Vp  and  Aejj  is  assumed  to  be 
linear,  i.e. 


r  —  Cpp)  (p2  Pi)  /  a  v 

I  max  —  K  rp  d/2  X  ^eff) 

where  Cps  and  Cpp  are  the  plate  suction  surface  and  pressure  surface  pressure  coefficients9. 
Aefj  is  set  equal  to  T  x  d/2  and  K  is  a  constant  of  proportionality  to  be  determined  from 
the  inviscid  CFD  results.  Once  K  is  fixed  the  value  of  Tmax  is  determined  purely  by  the 
cascade  geometry  and  density  wake  parameters. 

7.3  Solution  Procedure 

The  bound  vortex  circulation  strengths  Tk  at  each  discrete  vortex  panel  k  and  the  induced 
velocities  F *  from  the  forcing  terms  (described  shortly)  at  each  panel  3/4  chord  location 
are  related  by  an  aero-influence  coefficient  matrix  A  as  follows, 


{F}  =  [A){T} 

9{Cps  —  Cpp)  is  evaluated  at  the  plate  1/4  chord  location. 


(17) 
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The  A  matrix  stores  the  vertical  component  of  velocity  induced  at  the  3/4  chord  location 
((k)  of  each  panel  by  the  bound  vortex  of  strength  at  each  panel  1/4  chord  location  £(/;). 
The  A  matrix  can  be  obtained  by  substituting  appropriate  values  for  £  and  £  in  Equation  11 
for  a  given  cascade  space-chord  ratio  and  stagger  angle. 

The  forcing  vector  F *  contains  (1)  the  vertical  component  of  the  free  stream  velocity 
Un  and  (2)  the  velocity  induced  by  the  counterrotating  vortex  pairs  at  the  panel  3/4  chord 
locations.  The  F^  vector  is  constantly  updated  as  the  counterrotating  vortex  pairs  convect 
through  the  cascade.  The  velocity  induced  by  the  counterrotating  vortices  at  the  plate 
leading  edge  is  also  added  to  F^10. 

The  unknown  values  of  the  panel  bound  vortex  strengths  I\  can  now  be  evaluated  by 
inverting  the  aero-influence  coefficient  matrix  A  as  follows, 

{r}  =  [A]-1  {F}  (18) 


The  cascade  lift  (force  normal  to  the  flat  plate  along  the  —rj  direction)11  and  the  moment 
(clockwise  positive  about  the  plate  mid-chord)  is  then  obtained  using, 


N 

puY,v* 

(19) 

k=l 

N 

(20) 

k=i 


Note  the  influence  of  shed  vorticity  and  the  unsteady  contribution  to  lift  and  moment 
from  the  unsteady  circulation  dT^/dt  is  not  included  in  this  model.  The  remaining  sections 
describe  the  validation  of  the  model  and  the  parametric  studies  conducted  with  different 
density  wakes  and  cascade  geometries. 


7.4  Steady  State  Model  Validation 

The  steady  state  lift  coefficient  predicted  by  the  model  (i.e.  with  no  converting  density 
wakes)  is  shown  in  Figure  29  for  a  range  of  cascade  space-chord  ratio’s  and  stagger  angles. 
The  cascade  lift  coefficient  Clcascade  shown  here  is  normalized  by  the  lift  coefficient  of  an 
isolated  flat  plate  airfoil  Clairf0li.  A  single  vortex  panel  is  used  in  this  case.  All  curves  show 
an  asymptotic  value  of  1.0  as  the  space-chord  ratio  is  increased.  This  is  as  expected  and 
confirms  the  cascade  lift  coefficient  approaches  the  flat  plate  lift  coefficient  in  the  limit  of 

10 This  represents  the  influence  of  the  counterrotating  vortices  on  the  free  stream  flow, 

11  The  term  “lift”  will  be  used  subsequently  instead  of  “force”  to  help  indicate  this  direction. 
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Figure  29:  Cascade  model  results  for  the  change  in  cascade  interference  coefficient 
Clcascadt/Cla,rfo,i  with  space-chord  ratio  and  stagger  angle.  No.  of  panels  =  1. 


large  space-chord  ratios.  Also  notice  the  large  increase  in  cascade  lift  coefficient  for  stagger 
angles  /?i  >  45  degrees  and  space-chord  ratio  a  a  1.0.  This  is  due  to  the  increased  flow 
restriction  between  the  adjacent  flat  plates  at  these  high  values. 

The  corresponding  plot  obtained  using  a  conformal  mapping  method  by  Weinig  [9]  is 
shown  in  Figure  30.  The  overall  trends  compare  well  however  the  numerical  values  differ 
significantly  at  the  higher  stagger  angles.  Better  agreement  can  be  obtained  using  5  vor¬ 
tex  panels  to  represent  the  flat  plates.  The  corresponding  results  shown  in  Figure  31  now 
compare  to  within  2%  of  Weinig’s  results12.  The  model  results  for  the  steady  state  moment 
coefficient  about  the  plate  mid-chord  location  is  plotted  in  Figure  32  normalized  about 
the  moment  coefficient  for  an  isolated  airfoil.  An  equivalent  theoretical  prediction  could 
not  be  found  for  comparison  purposes  for  this  case.  The  moment  coefficient  nevertheless 
approaches  the  thin  airfoil  predicted  value  in  the  limit  of  large  space-chord  ratio. 

The  results  presented  here  confirm  the  model  can  adequately  predict  the  numerical 
values  and  trends  for  the  steady  state  cascade  lift  and  moment  coefficients.  The  model 
results  for  the  quasi-steady  lift  and  moment  profiles  induced  by  density  wakes  are  discussed 
next. 

12Results  for  stagger  angles  greater  than  60  degrees  could  not  be  obtained  using  this  model  since  the 
aero-influence  coefficient  matrix  .4  becomes  ill-conditioned  and  cannot  be  inverted. 
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Figure  30:  Weinig’s  conformal  mapping  prediction  for  the  cascade  interference  coefficient  k0 
{Cl cascade/ Cl atrf oil)  as  a  function  of  stagger  angle  7 eff  and  space-chord  ratio  (s/c)c//; 


Space-Chord  ratio  (T/c) 


Figure  31:  Cascade  model  results  for  the  change  in  cascade  interference  coefficient 
Cl  cascade /Cl  axrf oxi  with  space-chord  ratio  and  stagger  angle.  No.  of  panels  =  5. 


Figure  32:  Cascade  model  results  for  the  change  in  cascade  interference  coefficient 
CmcaScadtlCTnaxrfoii  with  space-chord  ratio  and  stagger  angle.  No.  of  panels  =  5. 


Property 

Value 

Space/chord  ratio 

Stagger  angle 

Max.  blade  thickness 

Blade  leading  edge  radius 

0.58 

15.0  deg. 

0.037c 

0.0075c 

Table  2:  Properties  of  the  NACA4F  cascade  geometry  used  for  the  inviscid  CFD  tests. 


7.5  Quasi-Steady  Model  Validation 

The  quasi-steady  model  validation  involves  comparison  of  the  model  prediction  for  the  den¬ 
sity  wake  induced  flat  plate  lift  and  moment  fluctuation  with  a  corresponding  inviscid  CFD 
prediction.  The  model  parameters  are  set  to  match  the  cascade  geometry  properties  and 
density  wake  properties  of  a  typical  inviscid  CFD  calculation.  For  the  cascade  properties,  a 
NACA4F  blade  geometry  was  chosen  with  properties  listed  in  Table  2.  For  the  density  wake 
a  width  of  0.1c  and  density  ratio  0.5  was  considered  a  representative  case.  The  model  pa¬ 
rameters  corresponding  to  this  set  of  cascade  and  density  wake  properties  is  listed  in  Table  3. 

The  model  results  and  the  inviscid  CFD  results  for  the  force  and  moment  fluctuations 
are  compared  in  Figure  33.  The  dashed  line  represents  the  model  prediction  and  the  solid 
line  represents  the  inviscid  CFD  prediction.  The  change  in  the  lift  and  moment  coefficients 
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Variable 

Value 

Number  of  vortex  panels  (N) 
Space/chord  ( s/c ) 

Cascade  stagger  angle  (/?i) 
Counterrotating  vortex  spacing  ( d ) 
Density  ratio  (P2/P1) 

Upstream  infl.  of  press,  field  (l\) 
Downstream  infl.  of  press,  field  (l?) 
Incidence  angle  (cr) 

Circulation  strength  constant  (K) 

20 

0.58 

15.0 

0.1c 

0.5 

0.3c  upstream  of  l.e. 

0.4c  downstream  of  l.e. 

5  deg.  relative  to  blade 
0.32 

Table  3:  Cascade  model  parameters  used  to  determine  the  flat  plate  force  and  moment  coefficient 
fluctuation  during  passage  of  a  density  wake  of  width  0.1c  and  density  ratio  0.5. 


are  shown  as  percentage  changes  from  the  steady  state  values.  Note  the  definition  for  the 
model  moment  coefficient  (positive  clockwise  about  the  mid-chord)  was  changed  to  posi¬ 
tive  counter-clockwise  about  the  mid-chord  for  comparison  purposes  with  the  inviscid  CFD 
results13.  The  maximum  fluctuation  in  the  model  lift  prediction  is  matched  to  the  inviscid 
results  by  adjusting  the  value  of  the  circulation  constant  K  (Equation  16).  20  vortex  panels 
are  used  to  discretize  the  cascade  flat  plates  for  this  case.  Use  of  more  than  20  panels  is 
found  to  change  the  amplitude  and  frequency  response  by  negligible  amounts. 

Significant  deviation  can  be  seen  between  the  model  and  the  inviscid  CFD  predictions. 
The  maximum  change  in  the  lift  response  occurs  at  a  time  lag  At  ss  0.4  after  the  corre¬ 
sponding  maximum  change  in  the  inviscid  CFD  lift  results.  Similarly  the  maximum  change 
in  the  moment  occurs  At  ~  0.2  before  the  corresponding  maximum  change  in  the  inviscid 
CFD  moment  results.  The  initial  decrease/increase  in  lift/moment  compares  well  how¬ 
ever  the  return  to  the  steady  state  values  does  not  compare.  In  addition  the  lift  peak  at 
t  =  0.1  predicted  by  the  model  has  a  larger  amplitude  and  occurs  later  than  the  similar 
peak  predicted  in  the  inviscid  CFD  results.  The  similar  peak  in  the  moment  response  is  not 
captured  by  the  model.  Furthermore  the  overall  maximum  change  in  the  moment  is  under 
predicted.  The  following  list  of  model  features  have  been  highlighted  as  possible  reasons  for 
the  observed  discrepancies. 

•  The  maximum  change  in  the  lift  and  moment  coefficients  occur  as  the  density  wake 
passes  over  the  blade  leading  edge  region.  The  blade  leading  edge  geometry  and 
corresponding  pressure  distribution  pattern  is  therefore  critical  to  the  overall  response 
profile.  The  model  however  replaces  the  finite  blade  radius  with  an  infinite  radius  and 
consequently  an  infinite  pressure  peak.  It  may  not  be  possible  therefore  to  adequately 
capture  the  initial  response  shape  with  a  flat  plate  geometry. 

3A1so  note  the  lift  coefficient  is  defined  normal  to  the  flat  plate  whereas  the  force  coefficient  obtained 
from  inviscid  CFD  simulations  is  defined  along  the  azimuthal  coordinate  direction.  The  difference  is  absorbed 
into  the  model  constant  K. 
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Figure  33:  Comparison  of  the  quasi-steady  model  results  and  the  inviscid  CFD  results  for  the 
fluctuation  in  the  blade  force  and  moment  coefficients  during  passage  of  a  density 
wake  of  width  0.1c  and  density  ratio  0.5.  A  Cl  =  (Cl  mar  Clm<an)/(Cl  mean))  AC*77i 
(Cm  max  Cmmean  )/(Cm  mean  )  • 


•  The  initial  rise  in  the  lift  response  curves  was  found  to  be  governed  by  the  induced 
velocity  field  of  the  clockwise  vortex  (vortex  B)  of  the  counterrotating  vortex  pair. 
The  clockwise  vortex  induces  an  upwash  velocity  at  the  plate  leading  edge  when  it 
is  located  downstream  of  the  plate  leading  edge.  This  results  in  an  increase  in  the 
flow  angle  of  attack  and  hence  an  increase  in  the  lift  coefficient.  The  model  moment 
response  does  not  capture  this  effect  however.  Since  the  moment  is  calculated  about 
the  plate  mid-chord  there  must  be  a  relative  difference  in  induced  velocity  at  the 
collocation  points  upstream  and  downstream  of  the  mid-chord  for  a  visible  fluctuation 
in  the  moment  to  occur.  The  induced  velocity  of  the  clockwise  vortex  at  the  leading 
edge  is  however  applied  equally  to  all  3/4  chord  collocation  locations  on  the  vortex 
panels.  There  is  little  if  no  net  difference  in  the  moment  therefore  due  to  the  increase 
in  flow  angle  of  attack.  The  induced  velocity  at  the  plate  leading  edge  should  therefore 
be  weighted  at  the  vortex  panel  collocation  locations  for  improved  comparison  with 
the  inviscid  CFD  results.  This  may  also  help  to  reduce  the  magnitude  of  the  initial 
increase  in  the  lift  coefficient. 

•  The  lag  between  the  model  and  inviscid  CFD  results  may  be  due  to  the  influence  of 
shed  vorticity  which  is  not  accounted  for  in  the  model. 


The  sensitivity  of  the  lift  and  moment  response  curves  to  changes  in  the  density  wake 
properties  and  the  cascade  geometry  is  examined  next. 
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7.6  Parametric  Results 

7.6.1  Lift  and  Moment  Sensitivity  to  Density  Wake  Properties 

The  NACA4F  cascade  geometry  is  used  as  the  baseline  geometry  for  the  sensitivity  anal¬ 
ysis.  The  circulation  constant  K  is  fixed  at  0.322.  This  sets  the  maximum  fluctuation  in 
lift  coefficient  to  7.3%  and  the  maximum  fluctuation  in  moment  coefficient  to  -6.8%  for  a 
density  wake  of  width  0.1c  and  density  ratio  0.5  (see  Figure  33). 

The  results  for  the  maximum  fluctuation  in  the  lift  and  moment  coefficients  with  vary¬ 
ing  density  wake  width  and  density  ratio  are  plotted  in  Figure  34.  The  variation  in  density 
ratio  is  represented  here  by  Marbles’  density  parameter  p*  (Equation  2).  The  model  data  is 
indicated  by  markers  joined  by  straight  lines.  The  corresponding  inviscid  CFD  results  are 
shown  superimposed  as  dashed  lines.  The  model  trends  agree  well  with  the  inviscid  CFD 
results.  In  particular  the  model  values  for  the  maximum  fluctuation  in  azimuthal  force  co¬ 
efficient  compare  to  within  2%  with  the  CFD  numerical  values  for  density  parameter  p*  <  0. 

It  is  interesting  to  note  how  the  model  predicts  the  non-linearity  in  the  maximum  lift 
and  moment  fluctuations  for  large  values  of  density  parameter  and  wake  width  similar  to 
the  non-linearity  observed  in  the  inviscid  CFD  results.  Initial  investigations  indicate  the 
non-linearity  to  originate  from  the  expressions  for  the  counterrotating  vortex  circulation 
strengths14  (Equation  14). 

7.6.2  Lift  and  Moment  Sensitivity  to  Cascade  Geometry 

The  density  wake  width  and  density  ratio  were  held  constant  at  0.1c  and  0.5  respectively 
as  the  cascade  stagger  and  space-chord  ratio  were  varied.  The  circulation  constant  K  was 
again  fixed  at  0.322. 

The  results  for  the  maximum  fluctuation  in  the  lift  and  moment  coefficients  with  varying 
cascade  space-chord  ratio  and  stagger  angle  are  plotted  in  Figure  35.  The  lift  coefficient 
is  more  sensitive  to  the  cascade  space-chord  ratio  than  to  the  stagger  angle.  This  may  be 
because  the  counterrotating  vortices  convect  at  the  same  stagger  angle  as  the  cascade15. 
The  reduction  in  lift  coefficient  fluctuation  as  the  space-chord  ratio  is  increased  is  most 
likely  due  to  the  increased  distance  between  the  counterrotating  vortices  and  the  blade  sur¬ 
faces.  The  moment  coefficient  fluctuation  also  decreases  with  increasing  space-chord  ratio 
and  shows  greater  sensitivity  to  the  cascade  stagger  angle. 


14  A  systematic  analysis  of  the  model  parameters  can  help  to  isolate  the  precise  source  of  this  non-linearity. 
15This  is  an  assumption  in  the  model 
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Figure  34:  Cascade  lift  and  moment  fluctuation  sensitivity  to  density  wake  width  w/c  and  density 
parameter  p*  (measure  of  density  ratio).  The  NACA4F  cascade  geometry  is  used  for 
all  tests.  Solid  lines  indicate  the  model  results.  Dashed  lines  indicate  the  inviscid 
CFD  results. 


Figure  35:  Cascade  lift  and  moment  sensitivity  to  cascade  stagger  angle  and  space-chord  ratio. 
Density  wake  of  width  0.1c  and  density  ratio  0.5  is  used  for  all  tests. 
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7.7  Conclusion 

A  model  to  predict  the  density  wake  induced  force  and  moment  fluctuations  on  cascade 
blade  rows  has  been  described  in  this  Section.  The  model  is  developed  using  the  following 
main  assumptions: 

1.  Effects  due  to  blade  camber  and  thickness  can  be  neglected  so  that  the  cascade  blades 
can  be  represented  as  flat  plates. 

2.  The  density  wake  can  be  represented  as  a  row  of  counterrotating  vortices  which  convect 
with  the  mean  flow. 

3.  The  influence  of  shed  vorticity  can  be  neglected. 

The  model  predicts  the  maximum  fluctuation  in  the  lift  and  moment  coefficients  to  (1) 
increase  with  density  wake  width,  (2)  increase  with  density  parameter  |p*|,  (3)  decrease  with 
increasing  space-chord  ratio  and  (4)  remain  relatively  insensitive  to  blade  stagger  angle. 

The  model  predicted  trends  for  the  force  and  moment  fluctuations  with  varying  density 
wake  properties  and  fixed  cascade  geometry  follow  the  same  trends  as  the  inviscid  CFD 
results.  Moreover  for  the  case  of  density  parameter  p*  <  0  the  numerical  values  for  the 
maximum  fluctuation  in  the  lift  coefficient  agree  to  within  2%  of  the  inviscid  CFD  predic¬ 
tions. 

The  model  results  for  the  moment  coefficient  however  under  predict  the  inviscid  CFD 
results  by  up  to  50%.  The  model  does  not  adequately  capture  the  extra  length  parameter 
(moment  arm)  required  to  determine  the  moment  coefficient  16 .  Thus  the  quasi-steady 
model  while  being  instructive  fails  to  match  the  unsteady  characteristics. 


16The  assumptions  in  the  model  regarding  the  magnitude  of  the  induced  velocity  at  the  vortex  panel 
collocation  points  needs  to  be  re-examined  to  determine  the  origin  of  this  difference. 
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8  Overall  Conclusions 

The  important  conclusions  from  this  study  are  listed  below. 

•  The  study  of  convecting  density  wakes  through  cascade  blade  rows  have  highlighted 
4  distinct  flow  mechanisms  which  generate  fluctuating  blade  forces  and  moments: 

1.  The  interaction  of  pressure  gradients  and  density  gradients  to  form  counterro¬ 
tating  vortex  pairs  and  concentrated  fluid  flux. 

2.  Axial  deflection  of  the  blade  passage  shock  wave  position. 

3.  Convecting  separation  bubbles. 

4.  Vortex  shedding. 

The  first  3  mechanisms  are  found  to  be  excited  by  the  passage  of  density  wakes,  the 
4th  mechanism  is  found  to  be  temporarily  suppressed  by  density  wakes. 

•  The  magnitude  of  the  blade  force  and  moment  fluctuations  generated  by  the  above 
mechanisms  are  found  to  be  characterized  by  a  set  of  non-dimensional  parameters 
based  on  density  wake  properties,  flow  Mach  number  and  blade  geometry: 

1.  Density  wake  width  w/c. 

2.  Density  parameter  p*. 

3.  Density  ratio  P2/P1. 

4.  Prandtl  Glauert  compressibility  factor  \/l  — 

5.  Blade  spacing  h/c. 

The  magnitude  of  the  blade  force  and  moment  fluctuations  have  also  been  quantified 
using  simple  functional  relationships. 

•  A  simple  cascade  flow  model  has  been  developed  to  establish  trends  for  the  density 
wake  -  cascade  blade  row  interaction.  The  model  captures  the  trends  for  the  variation 
of  the  blade  force  and  moment  coefficients  as  predicted  by  the  inviscid  CFD  flow  sim¬ 
ulations.  The  model  also  shows  the  maximum  fluctuation  in  the  force  and  moment 
coefficients  to  (1)  decrease  with  increasing  cascade  space-chord  ratio  and  (2)  remain 
relatively  insensitive  to  blade  stagger  angle. 
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14  Viscous  Flow  Solver 


A  two-dimensional,  time-accurate,  explicit  Navier-Stokes  flow  solver  developed  by  Hoy- 
ing  [11]  is  used  for  all  the  viscous  compressible  flow  simulations.  The  aim  of  this  section 
is  to  briefly  highlight  the  main  features  of  this  flow  solver  and  to  detail  the  computational 
grid  and  boundary  conditions  used  to  obtain  converged  solutions.  The  compressor  cascade 
properties  and  the  imposed  shape  of  the  density  wake  profiles  is  also  detailed.  The  char¬ 
acteristic  variables  used  to  non-dimensionalize  the  results  is  also  described.  Details  of  the 
inviscid  flow  solver  and  computational  grid  can  be  obtained  from  Ramer  [22]. 

14.1  Features  of  The  Viscous  Flow  Solver 

The  key  features  of  the  Navier-Stokes  flow  solver  are  summarized  below. 

1.  Finite  difference  discretization. 

2.  Fourth-order  spatial  accuracy  and  third-order  temporal  accuracy. 

3.  Dispersion  Relation  Preserving  (DRP)  scheme  (Tam  and  Web  [25]). 

4.  k-e  turbulence  model. 

5.  Wall  functions  (Chieng  and  Launder  [5]). 

6.  Non-reflecting  inlet  and  exit  boundary  conditions  (Giles  [7]). 

The  Dispersion  Relation  Preserving  scheme  used  in  this  solver  has  the  advantage  of  op¬ 
timizing  the  dissipation  and  dispersion  characteristics  of  the  finite  difference  method  over  a 
large  number  of  wavelengths  (than  say  a  four-stage  Runga-Kutta  scheme).  The  resulting 
high-order  scheme  also  reduces  the  computational  cost  by  requiring  fewer  grid  points.  As 
is  the  case  with  most  high-order  methods,  the  presence  of  high  frequency  waves  in  the 
solution  is  of  great  concern.  The  scheme  is  noted  as  being  able  to  capture  high  frequency 
waves  to  an  extent  that  unwanted  waves  will  remain  in  the  solution.  Artificial  damping 
based  on  the  fourth  derivative  of  pressure  at  each  node  however  was  found  to  be  effective 
in  removing  the  high-frequency  waves  without  creating  excessive  damping  of  the  desired 
waves.  Specific  implementation  details  are  described  in  Hoying  [11]. 

A  standard  two  equation  k-e  turbulence  model  was  implemented  in  the  solver  using  the 
following  model  constants, 

C v  =  0.09,  Ci  =  1.44,  C2  =  1.92,  <rK  =  1.0,  <r6  =  1.3  (21) 


These  constants  as  proposed  by  Launder  and  Spalding  [15]  are  appropriate  for  free  shear 
flows  and  flows  near  solid  boundaries.  A  two-stage  MacCormack  Predictor-Corrector 
method  is  chosen  to  integrate  the  k-e  model  equations.  The  flow  equations  and  the  k- 
e  model  equations  are  advanced  in  time  with  the  same  time-step. 
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Periodic  boundary 


Inlet  boundary 


7^ 

Exit  boundary 


Figure  36:  Schematic  of  the  computational  domain  and  boundary  conditions. 


Additional  computational  savings  are  achieved  by  using  wall  functions  to  model  the 
inner  portion  of  the  blade  boundary  layers.  Wall  functions  allow  the  first  near  wall  grid 
point  to  be  located  at  y+  «  30-150.  This  results  in  a  considerable  saving  of  the  number  of 
near  wall  grid  points.  The  wall  functions  calculate  the  wall  shear  stress  and  therefore  a  zero 
slip  velocity  at  the  wall  cannot  be  specified  as  a  boundary  condition.  A  finite  slip  velocity 
is  specified  instead  such  that  the  correct  wall  shear  stress  is  achieved  17 .  An  adiabatic  wall 
boundary  condition  is  also  applied  at  the  blade  surface. 

The  computational  domain  used  in  the  solver  is  shown  schematically  in  Figure  36.  The 
flow  upstream  and  downstream  of  the  blade  row  is  assumed  to  be  governed  by  the  Euler 
equations.  The  flow  in  the  blade  passages  surrounding  the  blades  is  assumed  to  be  gov¬ 
erned  by  the  Navier-Stokes  equations.  This  division  helps  to  reduce  the  computational  time 
since  the  solution  of  the  more  expensive  Navier-Stokes  equations  is  confined  to  a  physically 
smaller  region  of  the  flow. 


Non-reflecting  boundary  conditions  developed  by  Giles  [7]  are  implemented  at  the  inlet 
and  exit  boundaries  of  the  computational  domain.  The  aim  of  these  boundary  conditions 
are  to  allow  all  out-going  modes  of  the  solution  to  propagate  through  the  inlet  and  exit 
boundaries  without  reflecting  back  into  the  computational  domain.  Periodic  boundary  con¬ 
ditions  are  employed  along  the  side  boundaries  of  the  computational  domain  parallel  to  the 
axial  flow  direction.  These  boundary  conditions  are  shown  in  Figure  36. 


17 See  Hoying  for  more  details. 


Figure  37:  O-grid  detail  at  the  blade  leading  edge. 


14.2  Computational  Grid 

The  computational  domain  used  in  the  solver  consists  of  overlapping  H  and  O  grids  as 
shown  schematically  in  Figure  36.  The  H-grid  has  uniform  spacing  throughout  the  compu¬ 
tational  domain  with  a  cell  aspect  ratio  of  1.0.  The  orthogonality  of  the  H-grid  simplifies 
the  evaluation  of  numerical  derivatives  as  well  as  the  H  to  O  grid  interpolation  coefficients. 
The  O-grid  was  generated  using  an  elliptic  partial  differential  equation  method  developed 
by  Steger  and  Sorenson  [24].  The  O-grid  detail  near  the  blade  leading  edge  and  trailing 
edge  is  shown  in  Figure  37  and  Figure  38.  The  grid  cell  aspect  ratios  obtained  by  this 
method  are  close  to  1.0  near  the  blade  leading  and  trailing  edge.  This  is  a  typical  “rule  of 
thumb”  value  for  grid  cell  aspect  ratio  in  regions  where  the  flow  is  not  unidirectional  [20]. 
The  grid  is  also  clustered  closer  at  the  blade  surface  and  around  the  leading  and  trailing 
edge  to  help  resolve  the  higher  gradients  in  flow  quantities  in  these  regions.  The  Steger  and 
Sorenson  method  also  allows  grid  lines  at  the  surface  to  be  orthogonal  to  help  simplify  the 
evaluation  of  numerical  derivatives.  The  grid  is  near-orthogonal  in  regions  away  from  the 
surface. 

The  minimum  grid  resolution  required  to  adequately  capture  flow  features  is  a  primary 
concern  in  any  computational  study.  In  the  present  study  a  detailed  analysis  of  grid  de¬ 
pendency  was  not  conducted  however.  Instead  the  number  of  H-grid  and  O-grid  nodes 
was  increased  until  an  upper  limit  was  reached  for  the  time  required  to  achieve  a  converged 
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Figure  38:  O-grid  detail  at  the  blade  trailing  edge. 


solution  (i.e.  roughly  equal  to  3-4  weeks  on  a  DEC  Alpha  workstation).  The  resulting  0- 
grid  consists  of  421  nodes  around  the  blade  surface  and  32  nodes  radially  outward  from  the 
blade  surface.  The  H-grid  consists  of  161  nodes  in  the  azimuthal  direction  and  1283  nodes 
in  the  axial  direction. 


14.3  Cascade  Geometry  and  Blade  Profile 

The  blade  profile  used  in  this  study  is  the  General  Electric,  Low  Speed  Research  Com¬ 
pressor  (LSRC)  Stator-B  blade.  The  50%  span  section  of  this  blade  shown  in  Figure  39 
is  used  for  all  the  viscous  flow  simulations.  The  LSRC  blade  design  is  described  in  detail 
by  Wisler  [30].  The  cascade  geometry  used  for  the  viscous  CFD  simulations  are  listed  in 
Table  4. 


14.4  Density  Wake  Profile 

A  single  density  wake  was  convected  through  the  cascade  blade  row  in  each  simulation. 
The  density  wake,  density  variation  is  assumed  to  take  a  (1  —  cosine)  functional  form.  This 
ensures  a  smooth  transition  in  density  between  the  density  wake  and  the  free  stream  flow  to 
help  avoid  possible  numerical  instabilities  associated  with  step  changes  in  flow  properties. 
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Cascade  Geometry 

Value 

Chord/space  ratio 

Stagger  angle 

Max.  blade  thickness 

Blade  leading  edge  radius 

1.66 

24.5  deg. 

0.071c 

0.007c 

Table  4:  Properties  of  the  LSRC  cascade  geometry  used  for  the  viscous  CFD  simulations.  Ex¬ 
perimental  performance  results  from  Wisler  indicate  a  diffusion  factor  of  0.49  and  a  loss 
coefficient  of  0.01  for  this  blade  profile. 


Figure  39:  50%  span  section  of  the  LSRC  Stator-B  blade  used  for  the  viscous  simulations.  X 
denotes  the  location  about  which  moments  are  evaluated  (0.42c,  0.29c). 
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Figure  40:  Density  wake  profile,  w/c  =  0.2,  P2/P1  =  0.5. 


The  density  wake  profile  specified  to  represent  a  wake  width  0.2c  and  density  ratio  P2/P1  = 
0.5  is  illustrated  in  Figure  40. 

14.5  Converged  Solutions 

A  converged  baseline  solution  was  first  obtained  before  simulating  the  passage  of  density 
wakes  through  the  cascade.  Solution  convergence  was  determined  by  monitoring  the  time 
history  of  the  L2  norm  of  the  flow  variables  p,  pu,  pv  and  pe  where  p,  u,  v,  and  e  denote 
density,  axial  velocity,  azimuthal  velocity  and  energy.  A  typical  variation  in  the  L2  norm 
for  pe  is  plotted  in  Figure  41.  A  converged  solution  was  attained  when  the  L2  norm  for 
each  flow  variable  showed  no  further  significant  reduction  with  time.  In  addition  the  flow 
field  contour  plots  were  examined  to  ensure  that  unsteady  transients  were  no  longer  present 
in  the  solution. 


14.6  Non-Dimensionalization 

The  equations  of  motion  are  non-dimensionalized  using  the  following  reference  quantities. 

•  Length  :  blade  spacing  (s) 

•  Velocity  :  rotor  tip  speed  {Ur) 

•  Density  :  free  stream  density  (p0 Q) 

•  Temperature  :  free  stream  temperature  (Too) 
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Figure  41:  Convergence  history  for  the  L2  norm  of  density  p  x  energy  e.  Moo  —  0.53. 


•  Viscosity  :  free  stream  viscosity  {poo) 


•  Pressure  :  PooU\ 


•  Energy  :  U7R 


For  the  cascade  geometry  considered  in  this  study  the  blade  chord  c  and  inlet  total 
velocity  Uoo  are  more  relevant  non-dimensional  variables  for  the  length  and  velocity  scales18. 
These  variables  are  therefore  used  in  the  definitions  for  Reynolds  number,  Mach  number, 
force  coefficients,  moment  coefficient  (positive  clockwise  about  x/c  =  0.42,  y/c  —  0.29)  and 
skin  friction  coefficient  as  follows, 


18 


An  additional  factor  of  0.5  for  the  pressure  non-dimensionalization  is  also  more  appropriate. 
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Note  the  flow  Reynolds  number  and  Mach  number  based  on  the  inlet  total  velocity  and 
blade  chord  are  determined  after  a  converged  solution  is  obtained  and  cannot  be  specified 
a-priori.  The  Reynolds  numbers  and  Mach  numbers  required  as  input  conditions  to  the 
flow  solver  are  based  on  the  blade  spacing  and  rotor  tip  speed  instead19. 


19Since  there  is  no  rotor  in  this  study  the  values  for  Reynolds  number  and  Mach  number  specified  as  input 
conditions  are  essentially  unit  values. 
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15  Baseline  Flow  Results 

The  following  plots  are  included  in  this  Section: 

•  Time  averaged  blade  pressure  distributions. 

•  Azimuthal  force  coefficient  Cy ,  axial  force  coefficient  Cx  and  moment  coefficient  Cm 
fluctuations  in  the  baseline  response. 

•  Time  averaged  boundary  layer  properties:  non-dimensional  momentum  thickness  0/c, 
non-dimensional  displacement  thickness  8* /c  and  skin  friction  coefficient  Cf. 

Note  the  following: 

•  The  moment  coefficient  is  defined  positive  clockwise  and  is  calculated  about  the  blade 
coordinates  xjc  =  0.42  and  y/c  =  0.29  (see  Figure  39). 

•  The  boundary  layer  properties  are  plotted  vs.  the  surface  distance  to  help  increase 
resolution  at  the  blade  leading  and  trailing  edge  regions. 
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Figure  44:  Time  averaged  blade  pressure  distribution.  =  0.63. 


Figure  45:  Time  averaged  blade  pressure  distribution.  M <*>  =  0.87. 
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Figure  46:  The  force  and  moment  coefficient  fluctuation  in  the  baseline  solution  of  Run  1.  Cy, 
Cx  and  Cm  are  the  blade  azimuthal  force,  axial  force  and  moment  coefficients  re¬ 
spectively.  r  =  convective  time  scale.  Moo  =0.15. 
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Figure  47:  The  force  and  moment  coefficient  fluctuation  in  the  baseline  solution  of  Run  2.  Cy , 
Cx  and  Cm  are  the  blade  azimuthal  force,  axial  force  and  moment  coefficients  re¬ 
spectively.  r  =  convective  time  scale.  Moo  —  0.53. 
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Figure  48:  The  force  and  moment  coefficient  fluctuation  in  the  baseline  solution  of  Run  3.  Cy , 
Cx  and  Cm  are  the  blade  azimuthal  force,  axial  force  and  moment  coefficients  re¬ 
spectively.  r  =  convective  time  scale.  Moo  =  0.63. 


Figure  50:  Time  averaged  boundary  layer  properties.  Moo  =  0.15. 


Figure  51:  Time  averaged  boundary  layer  properties.  Moo  =  0.53. 
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Figure  52:  Time  averaged  boundary  layer  properties.  Moo  =  0.63. 


Figure  53:  Time  averaged  boundary  layer  properties.  Moo  =  0.87. 
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16  Primary  and  Secondary  Response  Results 

The  following  are  included  in  this  Section: 

•  Parametric  test  results  for  the  (1)  Azimuthal  force  coefficient  Cy,  (2)  axial  force 
coefficient  Cx  and  (3)  moment  coefficient  Cm  fluctuations.  The  complete  list  of 
parametric  tests  are  detailed  in  Table  5.  Consecutive  plots  are  ordered  first  according 
to  the  density  ratio,  secondly  by  the  Mach  number  and  thirdly  by  the  wake  width. 

•  Functional  relationships  to  scale  the  maximum  fluctuation  in  the  force  and  moment 
response. 

Note  the  moment  coefficient  is  defined  positive  clockwise  and  is  calculated  about  the  blade 
coordinates  x/c  =  0.42  and  y/c  =  0.29  (see  Figure  39). 


Parameter 

Run  1 

Run  2 

Run  3 

Run  4 

Moo 

0.15 

0.53 

0.63 

0.87 

w/c 

0.1, 0.2, 0.4. 1.0 

0.1,  0.2,  0.4 

P2/P1 

0.25,  0.50,  0.75,  2.00 

Table  5:  Parametric  test  variables,  w/c  =  non-dimensional  wake  width,  P2/P1  =  maximum 
density  inside  wake  /  free  stream  density. 
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Figure  102:  Fluctuation  in  blade  force  and  moment  coefficients,  ui/c  =  1.0,  pi/ p\  =  0.25, 

Moo  =0.15. 

16.1  Functional  Relationships  -  Primary  Response 

The  change  in  the  maximum  fluctuation  of  the  primary  response  force  and  moment  coeffi¬ 
cients  can  be  expressed  by  the  following  algebraic  relationships, 


% 


%Cy  max  Cy  mean 

- ^ - 

^ Vmean 

Climax  Cmmean 
Criflmean 

Cx-max  CXmean 


% 


CxTi 


n  x  p*2  +  n2P*  +  n3 
n  \p*2  +  n2P*  +  n3 
nlP*2  +  n2p*  +  n3 


(22) 

(23) 

(24) 


The  constants  Ill,  U2  ^nd  Il3  in  the  functional  relationships  are  tabulated  in  Table  6 
through  Table  9  for  each  free  stream  Mach  number  flow.  Note  the  functional  relationships 
obtained  here  represent  “best  fit”  curves  to  the  data  in  a  least-squares  sense.  This  results 
in  a  small  non-zero  value  for  the  II3  constant,  often  2  orders  of  magnitude  smaller  than  the 
values  for  II 1  and  n2,  which  can  be  ignored. 

The  n2  constant  for  all  force  and  moment  coefficients  is  seen  to  increase  by  a  factor  of 
between  1. 5-2.0  as  the  wake  width  doubles.  The  value  for  ITx  also  increases  with  increasing 
wake  width,  however  in  a  much  less  defined  manner  by  a  factor  of  between  2-20  as  the 
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Table  6:  Constants  in  the  functional  relationships  for  the  maximum  fluctuation  in  the  force 
moment  coefficients  (primary  response).  Moo  =  0.15. 


Coefficient 

Wake  Width  (w/c) 

IIi 

Cy 

0.1 

8.7 

28.4 

0.2 

23.3 

53.0 

-0.1 

0.4 

49.9 

87.5 

1.0 

Cm 

0.1 

44.5 

72.2 

-0.3 

0.2 

1.0 

0.4 

87.3 

157.5 

0.4 

Cx 

0.1 

18.5 

1.4 

33.8 

60.8 

1.4 

0.4 

54.8 

97.3 

2.0 

Table  7:  Constants  in  the  functional  relationships  for  the  maximum  fluctuation  in  the  force  and 
moment  coefficients  (primary  response).  Moo  —  0.53. 


Coefficient 

Wake  Width  (w/c) 

IIi 

n2 

n3 1 

Cy 

0.1 

9.0 

29.7 

I'll 

33.2 

54.3 

0.4 

58.9 

90.7 

Cm 

■ 

75.3 

1.7 

52.4 

111.7 

99.1 

178.5 

Cx 

0.1 

10.9 

35.0 

34.1 

62.0 

pH 

0.4 

59.2 

97.3 

m 

Table  8:  Constants  in  the  functional  relationships  for  the  maximum  fluctuation  in  the  force  and 
moment  coefficients  (primary  response).  Moo  =  0.63. 


Coefficient 

Wake  Width  (w/c) 

IIi 

n2 

n3 

Cy 

0.1 

-10.6 

40.2 

0.4 

0.2 

-3.2 

63.5 

-0.2 

0.4 

62.0 

118.9 

-0.5 

Cm 

0.1 

-5.2 

71.7 

-1.3 

0.2 

-9.6 

100.6 

0.2 

0.4 

-2.1 

131.9 

6.3 

Cx 

0.1 

-21.3 

36.2 

0.5 

0.2 

-14.7 

60.7 

0.4 

0.4 

26.9 

95.6 

0.2 

Table  9:  Constants  in  the  functional  relationships  for  the  maximum  fluctuation  in  the  force  and 
moment  coefficients  (primary  response).  Moo  —  0.87. 
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wake  width  doubles.  The  largest  increase  in  ITi  is  seen  between  wake  widths  of  0.2c  and 
0.4c.  The  functional  relations  for  =  0.87  have  a  relatively  smaller  magnitude  for  IIi 
compared  to  the  lower  Mach  number  functional  relations.  This  is  due  to  the  near-linear 
variation  in  the  the  shock  wave  motion  (and  hence  blade  static  pressure  fluctuations)  with 
increasing  \p*\  alluded  to  in  the  previous  section. 

The  functional  relations  should  be  used  with  caution.  In  particular  for  density  parame¬ 
ters  p*  «  0  the  II3  constant  may  become  dominant  and  the  deviation  from  the  data  can  be 
large. 


16.2  Functional  Relationships  -  Secondary  Response 

Functional  relationships  can  be  derived  to  help  quantify  the  trends  in  the  secondary  re- 


sponse.  The  change  in  the  amplitude  of  the  force  and  moment  coefficients  with  density 
ratio  P2/P1  can  be  expressed  by  the  following  algebraic  relationships, 

07  Cy second 
Cy baseline 

—  ^l(P2/pl)2  +  $2(P2/pl)  +  $3 

(25) 

07  C'm,second 

°Cmbaseiine 

=  $l(p2/Pl)2  +  $2(P2/Pl)  +  $3 

(26) 

07  Cxsecond 
C^baseline 

=  $l(p2/Pl)2  +  $2{P2/Pl)  +  $3 

(27) 

The  constants  $1,  $2  and  $3  in  the  functional  relationships  are  tabulated  in  Table  10 
through  Table  13  for  each  free  stream  Mach  number  flow.  Note  the  functional  relation¬ 
ships  obtained  here  represent  “best  fit”  curves  to  the  data  in  a  least-squares  sense  and  are 
valid  only  for  P2/P1  <  TO-  Insufficient  data  points  are  available  to  obtain  good  curve  fits 
for  P2/P1  >  TO. 
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Coefficient 

Wake  Width  (w/c) 

*1 

$2 

$3 

Cy 

0.1 

18.0 

-32.0 

15.1 

0.2 

5.9 

-18.1 

13.0 

0.4 

9.2 

-25.7 

17.3 

1.0 

8.2 

-20.9 

13.8 

Cm 

0.1 

19.0 

16.1 

0.2 

10.3 

-25.2 

15.7 

0.4 

10.9 

-28.4 

18.4 

1.0 

0.7 

11.1 

Cx 

0.1 

8.8 

-16.1 

8.2 

0.2 

3.2 

-9.2 

7.0 

0.4 

8.7 

-19.9 

12.1 

1.0 

39.5 

-66.1 

28.1 

Table  10:  Constants  in  the  functional  relationships  for  the  maximum  fluctuation  in  the  force  and 
moment  coefficients  (secondary  response).  Valid  for  P2/P1  <  1-0.  Moo  —  0.15. 


Coefficient 

Wake  Width  (w/c) 

*1 

$2 

4*3 

Cy 

0.1 

11.7 

-21.5 

10.8 

0.2 

0.9 

-8.1 

8.2 

0.4 

-5.8 

1.7 

5.1 

Cm 

0.1 

3.2 

-7.9 

5.7 

0.2 

-5.4 

3.4 

3.0 

0.4 

-5.7 

4.2 

2.5 

Cx 

0.1 

6.8 

-13.9 

8.1 

0.2 

5.0 

-12.3 

8.3 

0.4 

-4.0 

-0.5 

5.5 

Table  11:  Constants  in  the  functional  relationships  for  the  maximum  fluctuation  in  the  force  and 
moment  coefficients  (secondary  response).  Valid  for  p-i/pi  <  1.0.  Moo  =  0.53. 
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Coefficient 

Wake  Width  (w/c) 

*1 

$2 

$3 

Cy 

0.1 

6.6 

-12.2 

6.5 

0.2 

3.2 

-9.8 

7.6 

0.4 

-3.1 

-1.0 

5.0 

Cm 

o.i 

3.9 

-8.8 

5.8 

0.2 

-0.2 

-5.0 

6.2 

0.4 

-1.5  ' 

-3.1 

5.5 

Cx 

0.1 

8.4 

-15.3 

7.9 

0.2 

4.6 

-11.7 

8.1 

0.4 

-1.9 

-4.1 

6.9 

Table  12:  Constants  in  the  functional  relationships  for  the  maximum  fluctuation  in  the  force  and 
moment  coefficients  (secondary  response).  Valid  for  p2/pi  <  1.0.  Moo  =  0.63. 


Coefficient 

Wake  Width  (w/c) 

*1 

$2 

$3 

Cy 

0.1 

-4.3 

2.2 

3.1 

0.2 

-6.4 

5.7 

1.7 

0.4 

4.2 

-9.3 

6.2 

Cm 

0.1 

-2.5 

0.7 

2.8 

0.2 

-5.5 

5.0 

1.4 

0.4 

-1.2 

-0.8 

3.0 

Cx 

0.1 

-0.9 

-1.1 

3.0 

0.2 

-3.9 

3.3 

1.7 

0.4 

6.2 

-11.3 

6.2 

Table  13:  Constants  in  the  functional  relationships  for  the  maximum  fluctuation  in  the  force  and 
moment  coefficients  (secondary  response).  Valid  for  P2/ pi  <1.0.  A/00  =  0.87. 
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17  Compressibility  Scaling  Of  The  Maximum  Force  and  Mo¬ 
ment  Fluctuations 

This  Section  includes  the  following  plots: 

•  The  maximum  fluctuation  in  the  (1)  azimuthal  force  coefficient  Cy ,  (2)  axial  force  co¬ 
efficient  Cx  and  (3)  moment  coefficient  Cm  (primary  response)  scaled  by  the  Prandtl- 
Glauert  compressibility  factor  \/l  —  M^.. 

Note  the  results  for  the  Mx>  =  0.15  calculation  are  overlayed  in  each  plot  for  comparison 
purposes. 
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Figure  108:  Comparison  of  the  Moo  =  0.15  and  the  Moo  = 
fluctuation  in  the  azimuthal  force  coefficient. 
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Figure  109:  Comparison  of  the  Moo  =  0.15  and  the  Moo  =  0.87  viscous  results  for  the  maximum 
fluctuation  in  the  azimuthal  force  coefficient. 
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Figure  114:  Comparison  of  the  A/oo  =  0  15  viscous  results  and  the  inviscid  results  for  the  maxi¬ 
mum  fluctuation  in  the  moment  coefficient. 
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Figure  115:  Comparison  of  the  A/oc  =  0  15  and  the  A/oo  =  0.53  viscous  results  for  the  maximum 
fluctuation  in  the  moment  coefficient. 


Figure  116:  Comparison  of  the  Moo  =  0.15  and  the  Moo  =  0.63  viscous  results  for  the  maximum 
fluctuation  in  the  moment  coefficient. 


Figure  117:  Comparison  of  the  Moo  —  0.15  and  the  M0 o  =  0.87  viscous  results  for  the  maximum 
fluctuation  in  the  moment  coefficient. 


s 


14 


18  Cascade  Model  :  Matlab  Source  Code 


Filename:  cascade. m 


%  [Cascade  flow  model  v7.0] 

*/,  Written  by  Sanith  Wijesinghe 
7.  MIT  Gas  Turbine  Lab.  July  1998. 

#/. 

'/.  Matlab  script  file. 

y. 

*/,  Lift  and  moment  fluctuation  on  a  flat  plate  cascade  induced  by  a 
*/,  convecting  density  wake. 

y. 

’/,  Notes: 

y,  o  Cascade  flat  plate  extends  from  -0.5  to  +0.5. 

*/,  o  Uncomment  lines  where  indicated  for  either  steady  or  quasi-steady 
*/,  simulation. 

*/.  List  of  Symbols: 


%  N 

y.  z 
y.  x 

*/,  alpha 
*/,  betal 
*/.  sig 


=  no.  of  vortex  panels 
=  vector  of  bound  vortex  locations 
=  vector  of  collocation  locations 
=  flow  angle  of  attack 
=  stagger  angle  of  blades 
=  gap-chord  ratio 


’/,  beta2 
y.  gamax, 

y.  d 

•/•  r2rl 
•/•  11 
y.  12 

y.  sa 

y.  sb 

X  k 


=  stagger  angle  of  counterrotating  vortices 
gbmax  =  max.  circulation  of  the  counterrotating  vortices 
=  separation  of  counterrotating  vortices 
=  density  ratio  rho_2  /  rho_l 
=  max.  upstream  extent  of  press,  field  } 

=  downstream  extend  of  press,  field  >  defines  tanh  func. 
-  location  of  vortex  A  (counterrotating  pair) 

=  location  of  vortex  B  (counterrotating  pair) 

=  empirical  constant  determined  from  inviscid  CFD 


y.  cl  =  lift  coeff.  normal  to  flat  plate  {L/(0.5*rho*U~2*c)} 

*/.  cm  =  moment  coeff.  about  mid  chord  {M/(0. 5*rho*U*2*c"2) } 

y. 

y. - 


clear  all; 
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*/,  Initialize 


global  xl 

sig  betal 

7.7.  Set  global  parameters 

alpha 

= 

5*pi/180; 

7.7.  General  parameters 

P 

= 

0.0; 

q 

= 

0.0; 

veil 

= 

0.0; 

vel2 

= 

0.0; 

flagl 

= 

0.0; 

N 

= 

5;  7.7. 

Cascade  parameters 

sig 

= 

0.666; 

betal 

= 

20*pi/180; 

t 

= 

1 :  N ; 

z 

= 

-0.5  +  (t-1) 

.*(1/N)  +  (1/N)*(0.25) ; 

X 

= 

-0.5  +  (t-1) 

. *(1/N)  +  (1/N)*(0.75) ; 

11 

= 

-0.80; 

7.7.  Wake  parameters 

12 

= 

-0.10; 

n 

= 

-2; 

beta2 

= 

betal ; 

d 

= 

0.1; 

r2rl 

= 

0.50; 

k 

= 

0.322; 

coef f a 

= 

6/(12-11); 

coef fb 

= 

-3* (11+12)/ (12-11) ; 

7,  Start  Main  loop 

7,++  Uncomment  this  section  for  steady  flow  simulation 

7,  for  betal  =  0*pi/180:10*pi/180:60*pi/180  7.7.  Step  in  stagger  angles 

7.  p=p+l 

7,  for  sig  =  0.01:0.1:4.0;  7.7.  Step  in  gap-chord  ratios 

7.  q=q+l ; 


7.++ 

7.  Evaluate  aero-influence  coefficient  matrix 
for  a=l :N, 
for  b=l :N, 

xl  =  x(a);  7.7.  Set  collo.  pt.  in  kernel  function, 

suml  =  -imag(kernel(z(b))) ;  7.7.  B.  vortex  z(b)  infl.  at  loc.  x(a)  . 
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A(a,b)  =  suml; 


7.7.  A  matrix 


end 

end 


7.++  Uncomment  this  section  for  quasi-steady  flow  simulation 

7.  Convect  vortex  pair 
P  =  P  + 1 

for  n  =  -1.0:0.02:2.5; 
q  =  q  +  1; 


7.++ 


sa  =  n; 
sb  =  sa  -  d; 

7.  Hyperbolic  tangent  variation  in  vortex  strength 

if  sa  <=  11 
ga  =  0.0; 

elseif  sa  >  11  &  sa  <=  12 

ga  =  gamax*(l+tanh(coeffa*sa  +  coef fb) ) *0 . 5 ; 
elseif  sa  >  12 
ga  =  gamax; 
end 

gb  =  -ga; 

7.  RHS  source  terms 

for  v=l:N  7.7.  Counterrotating  Vortex  pair  infl.  at  panels  +  l.e. 

pa  =  imag(i*ga*exp(i*beta2)*coth((pi*(l/sig)*(x(v)-sa) 

+  i*pi/2)*exp(i*beta2))) ; 

pb  =  imag(i*gb*exp(i*beta2)*coth((pi*(l/sig)*(x(v)-sb) 

+  i*pi/2) *exp(i*beta2) ) ) ; 

la  =  imag(i*ga*exp(i*beta2)*coth((pi*(l/sig)*(-0.5-sa) 

+  i*pi/2)*exp(i*beta2))) ; 

lb  =  imag(i*gb*exp(i*beta2)*coth((pi*(l/sig)*(-0.5-sb) 

+  i*pi/2) *exp(i*beta2) ) ) ; 
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w(v)  =  pa+pb+la+lb-s in ( alpha) ; 
end 


gam 

=  inv(A)*w’;  ! 

7,7,  Solve  for  Gamma 

cl 

=  4*sig*(sum(gam)) ;  ! 

7,7,  Calculate  lift  coeff. 

cm 

=  4*sig*(gam,*z’) ;  ! 

7,7,  Calculate  moment  coeff. 

kint 

=  cl/(-2*pi*sin(alpha)) ;  ! 

7.7.  Interf.  ->  cl(cas . ) /cl (air . ) 

*/,  Solution  vector 


outl(p,q) 
out2(p,q) 
out3(p ,q) 
out4(p,q) 
out5(p,q) 
out6(p,q) 


=  n; 

=  cl; 

=  cm; 

=  kint; 

=  sig; 

=  betal; 


7.  Calculate  gamax  and  gbmax 


if  flagl  ==0.0 
7,  gamma (x)dx  ->  gamma (x) 
gl  =  gam/(l/N) ; 

7,  Set  gamax  and  gbmax  using  gamma<30.25  =  (Cps-Cpp)  and  r2/rl 

gamax  =  -k* (gl (5) )*(r2rl-l .0) 
gbmax  =  k*(gl(5))*(r2rl-l .0) 

flagl  =  1.0; 
end 


7.++  Uncomment  this  section  for  steady  flow  simulation 

7.  end 
7.  q=0; 


7.++ 

end 

7.  Maximum  change  in  the  lift  and  moment  coefficients 
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chcl  =  -((max(out2(l, :))-out2(l,l))/out2(l,l))*100 
chcm  =  ( (max(out3(l , : ) )-out3 (1 , 1) )/out3(l , 1) ) *100 

7, - <END> - 


Filename:  kernel. m 


7.  Kernel  function  to  be  used  with  cascade  model. 

% 


% 

y. 

y. 


Written  by  Sanith  Wijesing 
MIT  Gas  Turbine  Lab.  July  1998. 


function  [t]  =  kernel (e) 

global  xl  sig  betal 

nu  =  size (e , 2) ; 
sx  =  xl*ones(l,nu) ; 
ztl  =  sx-e; 

t  =  (-i*exp(i*betal) *coth(pi* (ztl) * (1/s ig) *exp(i*betal) ) ) ; 


<END> 
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